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ABSTRACT 

We present maps, source catalogue and number counts of the largest, most complete 
and unbiased extragalactic submillimetre survey: the 850 fim SCUBA HAlf Degree 
Extragalactic Survey (SHADES). Using the Submillimetre Common-User Bolometer 
Array (SCUBA) on the James Clerk Maxwell Telescope (JCMT), SHADES mapped 
two separate regions of sky: the Suha,Tu/ XMM-Newton Deep Field (SXDF) and the 
Lockman Hole East (LH). Encompassing 93 per cent of the overall acquired data 
(i.e. data taken up to 2004 February 1), these SCUBA maps cover 720arcmin^ with 
an RMS noise level of about 2 mJy and have uncovered > 100 submillimetre galaxies. 
In order to ensure the utmost robustness of the resulting source catalogue, data re- 
duction was independently carried out by four sub-groups within the SHADES team, 
providing an unprecedented degree of reliability with respect to other SCUBA cata- 
logues available from the literature. Individual source lists from the four groups were 
combined to produce a robust 120-object SHADES catalogue; an invaluable resource 
for follow-up campaigns aiming to study the properties of a complete and consistent 
sample of submillimetre galaxies. For the first time, we present deboosted flux densi- 
ties for each submillimetre galaxy found in a large survey. Extensive simulations and 
tests were performed separately by each group in order to confirm the robustness of 
the source candidates and to evaluate the effects of false detections, completeness, and 
flux density boosting. Corrections for these effects were then applied to the data to 
derive the submillimetre galaxy source counts. SHADES has a high enough number 
of detected sources that meaningful differential counts can be estimated, unlike most 
submillimetre surveys which have to consider integral counts. We present differential 
and integral source number counts and find that the differential counts are better fit 
with a broken power-law or a Schechter function than with a single power-law; the 
SHADES data alone significantly show that a break is required at several mJy, al- 
though the precise position of the break is not well constrained. We also find that an 
850 fim survey complete down to 2 mJy would resolve 20-30 per cent of the Far-IR 
background into point sources. 

Key words: submillimetre - surveys - cosmology: observations - galaxies: evolution 
- galaxies: formation - galaxies: starburst - galaxies: high-redshift - infrared: galaxies 



1 INTRODUCTION 

Deep blank-field surveys conducted with the Submil- 
limetre Common-User Bolometer Array (SCUBA) on 
the James Clerk Maxwell Telescope (JCMT) have 
resolved as much as 50 per cent (depending on 
the survey depth) of the Far-Infrared Background 
(FIB) into discrete, high redshift-s ources with flux den- 
sity levels of 5*85 ^ 2 mJy iSraail. lyison fc Blainl 

1997 , iHuehes et aT 'l998', 'Barecr ct al.' '1998*. 'Blain et_alj 
19991 . iBarger. Cowie fc Sand ers 1999. Bales et al. 20001 



Cowie. Barger fc KneihibooilScott et al.il200llWebb et alJ 



20031. iBorvs et alJl2003ll . Intensive campaigns to study these 



submillimetre galaxies (SMGs) and millimetre galaxies at 
other wavelengths have shown that they are primarily pow- 
ered by star- formation, a lthough many do harbo ur an Active 
Galactic Nucleus (AGN: [Alexander et al]l2005^ . 

The implied star-formation rates are extremely high 
(100-1000 M0yr-i), and since the well-known local 



JCondonI Il992 


) apparently 


sources (e.g 


Kovacs et alJ 



120061) . the high submillimetre luminosities result in cor- 
respondingly high radio luminosities, which are detected 
in deep 1.4 GHz images. The resolution afforded by radio- 
interferometers results in precise optical identiflcations of 
^ 50 percent o f all known SMGs, and in this manner 
IChapman et alJ (12005,1 derive a median redshift of (z) ~ 2.2 



for the population using deep Keck spectroscopy. While 
much of what is kno wn about SMGs is based on the 
radio detected subset ilvison et al.l l2002l) . other studies 
(e.g. Pope ct al. 200^ find no significant differences in the 

radio-undetected population. 

Many investig a tions (e.g. [S niail^ Ivison fc Blainl 



Barger et al.l 12001 llvison et alj 120021: ISmail et alJ 



1997 



2002 



Webb et al."200^:'Borvs ct al' 2004: 'W ang. Cowie fc Barged 
2004:.Greve et al. 20 04) have suggested that SMGs are Ukely 



to be associated with an early-phase in the formation of 
massive galaxies. The intensity of their starbursts, the re- 
sulting high metallicity, along with their large dynamical 



masses, high 


gas 


fractions and inferred 


strong clustering 


ifCreve et al. 


2005; Swinbank et al. 


2004; Blain et alJl2004^ 



arc all suggestive of a close link to the formation phase of 
the most massive spheroids (e.g. ISmail et al.ll2004ll . 

Despite this progress, the samples being used in any 
given study are typically quite small, since SCUBA could 
only detect roughly one SMG per night in good weather 
conditions. And while the total number of SMGs detected 
over the lifetime of the instrument^ is now ~ 400, these are 
mainly drawn from small {<, 100 square arcminute) fields 
spread all over the sky, each observed and reduced by dif- 
ferent groups using different techniques and source identi- 



SCUBA was operational between 1997 and 2005 
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fication criteria. The desire to obtain a well characterised 
sample of hundreds of SMGs in a large, contiguous area was 
the motivation for the SCUBA HAlf-Degree Extragalactic 
Surve y (SHADES: iMortier et"al]l2005l Ivan Kampen et alJ 
I2OO5I) . 

SHADES is an ambitious wide extragalactic submil- 
limetre survey, split evenly between 2 separate regions of 
sky: the Suha.ru/ XMM-Newton Deep Field (SXDF) and the 
Lockman Hole East (LH). The aim was to map 0.5 deg^ (~ 7 
times the area of its predecessor, the SCUBA 8-mJy Survey; 
IScott et al.ll2002^ to a comparable depth of ct ~ 2mJy at 
850 ^m, which is roughly three times the 1 a confusion limit 
imposed by the underlying s ea of fainter unresolved sources 
in the coarse SCUBA be am (Hose 2001: H ughes et alJll998l: 
ICowie. Bareer fc KneibI 120021 . The survey began in late 
2002, and finished when SCUBA was decomissioned in late 
2005. SHADES is the largest survey in terms of observing 
time carried out with SCUBA. 

Many scientifically powerful results from SHADES will 
come from comparisons of the properties of sources found 
in the SHADES catalogue with observations at other wave- 
lengths. Some of these include radio identifications (Ivison 
et al. in preparation), FIR- radio photometric redshift esti- 
mates (Aretxaga et al. in preparation), and submillimetre- 
Spitzer-based SEDs and photometric redshifts (Clements et 
al., Eales et al. and Serjeant et al. in preparation). In this 
paper we present new constraints on the numbers of submil- 
limetre sources as a function of 850 fim flux density based on 
the SHADES catalogues. A complementary analysis effort to 
constrain the numbers of faint submillimetre sources directly 
from our maps without the intermediate st ep of making a 
catalogue, a so-called P{D) approach (e.g. ICondonlll974l 
IScheueJlToTl . will be reported separately. 

The main results presented in this paper are the 
SHADES catalogue and 850 /im number counts. The ob- 
servations are summarised in Section |5| Section |3 provides 
details about the data reduction. For the first time in a 
submillimetre survey the data are processed by four inde- 
pendent data reduction pipelines in order to increase the 
robustness of the results. In Section0J the amalgamation of 
the four source lists into a master 850 /im SHADES cata- 
logue is described. In Section |^ two different approaches to 
derive the number counts are described, compared and con- 
trasted. The SHADES difi'erential count measurements are 
provided here. In Section |H| models are fit to the differen- 
tial counts, and the cumulative source counts are computed 
in order to compare them with previous data. Concluding 
remarks regarding what was learned about data analysis by 
comparing the different reductions, as well as advice for fu- 
ture surveys, are given in Section |7| 



2 OBSERVATIONS 

The Suhsiiu/ XMM-Newton Deep Field (SXDF) and Lock- 
man Hole East (LH) regions are centred at (J2000) 
RA=2''17™5r:5, Dec=-5°00'18'.'5 and RA=10''52'"26':7, 
Dec=57°24'l2'.'6, respectively. They were observed with a 
resolution of 14.8 arcsec an d 7.5 arcsec at 850 a nd 450 fim, 
respectively, with SCUBA (iHoUand et al.lll999l) on the 15- 
m JCMT atop Mauna Kea in Hawaii. A detailed account 
of the survey design and observing strategy is given in 



iMortier et alJ (l2005ll but is summarised briefly here. The 
observing strategy consists of making three overlapping jig- 
gle maps for each of six different chop throws (30, 44 and 
68 arcsec) and chop position angle (PA, and 90 deg in 
RA/D ec coordinates) combinations, motivated bv lEmersonI 
il995h . Each set of six observations spans a range of air- 
mass values, while the observations are limited to a 225 GHz 
atmospheric opacity range of 0.05 < rcso < 0.1, so as 
to maintain uniform coverage over the entire survey area. 
The SHADES survey data collection took place between 
1998 March and 2005 June (including the epoch of exist- 
ing data fr om the LH portion of the SCUBA 8-mJy Survey, 
IScott et al . 2002, Fox ct al. 2002). 

After a series of technical problems, SCUBA was of- 
ficially decommissioned in 2005 July, with only a modest 
amount of SHADES data taken in its last year. Only SCUBA 
data taken up until 2004 February 1, covering a total area of 
720arcmin^ down to an RMS noise level of ~ 2.2 mjy, are 
included here. The sky coverage is approximately evenly dis- 
tributed between the two SHADES fields. Complete SCUBA 
SHADES maps will contain only an additional roughly 7 
per cent of the planned area, although these will be com- 
plemented with 1.1mm AzTEC data ( Wilson et al...20oi) 
which will be presented in a future paper. 



3 DATA REDUCTION: FOUR INDEPENDENT 
SHADES PIPELINES 

In order to ensure the utmost robustness of the result- 
ing 850 fim source catalogue, data reduction was indepen- 
dently carried out by four sub-groups drawn from within the 
SHADES consortium; we refer to these as Reductions A, B, 
C and D. This provides an unprecedented degree of reliabil- 
ity with respect to other SCUBA source catalogues available 
from the literature. The reduction of the 450 data, which 
are only of limited use, is discussed in Appendix IB3I A ma- 
jor strength of the SHADES analysis comes from the con- 
cordance of these reductions which have not been modified 
to bring them into agreement. 

All groups perform the basic reduction steps of com- 
bining the nods (i.e. re-switching), flat-fielding, extinction 
correcting, despiking, removing the sky signal, calibrating, 
rebinning the data into pixel grids, and finally extracting 
sources. All reduction groups also correct the astrometry 
of the SXDF observations prior to December 26, 2002 by 
—6.7 arcsec in Declination, in order to account for a JCMT 
pointing catalogue error for oCeti. For all data taken on 
or after 2002 December 6 the known dead bolometer 'G9' 
was manually removed from further analysis. Groups dealt 
with the presence of a 16 sample noise spike in the data 
in different ways , as fir st noted by Bor ys ct al. ( 2004 ) and 
ISawicki fc WebbI ^2005^ (see also iMortier et alJl2005r : this 
is further discussed in Appendix 1X1 



3.1 Description of the four independent 
reductions 

Reduction steps were performed by each group indepen- 
dently using some com bination of SURF (S CUBA User Re- 
duction Facility; .Jenness fc Lightfootlll993) and their own 
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locally-developed code s. The basic steps of the analysis pro- 
cedure are described in lMortier et alj i2005l) . Each reduction 
group has made different choices with respect to various de- 
tailed steps, which are described in Tabled Several of the 
choices merit attention and are discussed below. 

The teams each produce signal-to-noise ratio (S/N) 
maps convolved with the telescope point spread function. 
This procedure is mathematically equivalent to fitting the 
point spread function (PSF) to the data centred on each 
pixel and is the minimum variance estimator of the bright- 
ness of i solated point source s if the background noise is white 
(see e.g. Ivon Hoerneilll96'iD . In all four reductions, flux den- 
sity measured when bolometers in the 'off-beam chop' of the 
telescope point at the region in question has essentially been 
folded in to the likelihood maps in order to increase the sen- 
sitivity of the maps. Two approaches have been used here: 

(1) fitting each pixel in the map to the multi-beam PSF; or 

(2) folding in the flux from the off-beams in the timestream 
data and then fitting each pixel in the map to a single Gaus- 
sian. After calibrating, the final map noise RMS values for 
the LH and SXDF fields are 2.1 and 2.0mJy beam"^, re- 
spectively. 

Finally, the presence of the 16 sample noise spike does 
not warrant any major pipeline alterations, since its effect on 
these data was found to be almost negligible; an account of 
the investigation into this issue is presented in Appendix IXI 
Several different checks were performed to test the Gaussian- 
ity of the noise and to check individual source robustness by 
examining spatial and temporal fits to point sources in the 
maps and this analysis is presented in Appendix lEl 



(monthly average FCF values have been tabulated on the 
JCMT website^, and were used to calibrate the data in Re- 
duction D). Over-reliance on individual FCF measurements 
could inject noise into the analysis. Ifowever, Reduction B, 
which makes the most use of measured FCF values, produces 
maps whose noise is as low as the other reductions. 

3.4 Pixel size 

Reductions A and B make 'zero-footprint' maps using an 
optimal noise-weighted drizzli ng algorithm (essentiall y the 
limiting form of the metho d of|Fru chtcr fc Ifook 2002") with 
a pixel size of 1 arcsec^ (see lMortier et al . 2005). Reductions 
C and D use 3x3 arcsec pixels. The p ressure to use larger 
pixels is driven by the observation of iBorvs et alJ ll2003fl 
that statistical analysis becomes unpredictable when there 
are too few observations per pixel. The concern that large 
pixels might lead to larger uncertainty in the positions of the 
detected sources is not supported in the data (as discussed 
in Section HTfl . 



3.2 Opacity 

The dominant cause of variation in the multiplicative factor 
between detector voltage and inferred source flux density is 
temporal variation in the atmospheric opacity. In all four 
reductions detector voltages are div ided by the atmospheri c 
opacity inferred, using coefficients in lArchibald et all ^2002^ . 
from either the 225 GHz rcso measureme nts or the 18 3 GHz 
line of sight water vapour monitor fWVM: IWiedneJl 998l at 
the JCMT. The differences in strategy described in Tabled 
do not lead to measurable differences in the data. Largely 
because of SHADES, there is a substantial body of data cor- 
relating opacity measured via secant scans ('sky dips') for 
various chop-throws with the opacity inferred from measure- 
ments of the rcso or the JCMT WVM. Opacity monitors 
provide a lower noise measurement of opacity than is pro- 
vided by any single sky-dip and do not cost any observing 
time (and are therefore preferred). 



3.3 Calibration 

The flux conversion factor (FCF) converts opacity-corrected 
bolometer voltages to flux densities. It is essentially the re- 
ceiver responsivity. The FCF is fairly stable on both short 
and long time scales, but might depart from the long term 
average as the telescope mirror is distorted by thermal gra- 
dients near to sunset and sunrise (Jenness et al...2002l'l . Be- 
cause the FCF at 850 /^m is so stable, the long term average 
might be a better estimator of the FCF at a given moment 
than any single measurement taken during stable conditions 



http: / / www.jach.hawaii.edu/JCMT / continuum / calibration / sens / gains.html 
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Table 1. Data Reduction Procedures. 'Secondary Extinction Monitor' refers to the source of the estimation of the sky opacity when the 
'Primary Extinction Monitor' was unavailable. 



Step 


Reduction A 


Reduction B 


Reduction C 


Reduction D 


Reduction 
Code Used 


IDL-based routines used in 
the SCUBA 8-mJy Sur- 
vey IScott ct al. 2003). See 


Reduction A's IDL pipeline 
altered at the extinction 
correction, calibration and 


SURF scripts for reswitch- 
ing and fiat-fielding and lo- 
cally developed code from 


SURF scripts for roswitch- 
ing and fiatfielding and 
locally developed C code 




ISerieant et al.« Ii2003r. 


source extraction phases. 


IChaDinl J200J) thereafter. 


thereafter jBorvs et alJ 






SeelMortier et alj J2005D. 




120031). 


Primary 

Extinction 

Monitor 


Polynomial fit to rcso- 


Interpolated WVM read- 
ings. 


Same as Reduction A. 


Nearest WVM reading. 


Secondary 
Extinction 
Monitor 


Linear fit to neighbouring 
sky-dips and interpolation. 


TCSOi sky-dip values. 


Interpolated sky-dips. 


Nearest rcso reading, lin- 
ear fit to neighbouring sky- 
dips and interpolation. 


Cosmic 
Rays 


Iterative 3 cr cuts until no 
signal is removed. 


~ 32 successive 3 cr cuts. 


One 3 cr cut. 


Iterative 4.5 (T cuts until no 
signal is removed. 



Baseline 


Subtract the temporal 


Subtract mode signal iter- Subtract array median. 


Subtract mean of non- 


Subtrac- 


modal sky level obtained 


atively with cosmic rays. 


noisy bolometers from all 


tion 


from a fit to all bolometers 


Handle bolometers exhibit- 


bolometers iteratively with 




in the array iteratively 


ing excess 16 sample noise 


cosmic ray removal. 




with cosmic ray removal. 


separately. 





Calibration Calculate an FCF for each 
half night and apply to the 
data. 



Linearly interpolate be- 
tween all stable FCF mea- 
surements during a night, 
selected by eye. 



Mean of measured FCF 
before and after sun- 
rise/set (i.e. two averages 
per shift). Each of 3 chop 
throw amplitudes handled 
separately. 



Monthly average 
FCFs (single value of 
219 JyV-^beam-l for 
data before January 1999). 



Flux Den- Inverse variance weighted 
sity Maps fiux density summed in 
1 arcsec pixels. Each chop 
throw mapped separately 
(six flux density maps) and 
then combined to produce 
a single maximum likeli- 
hood flux density map. 



Same as Reduction A. 



Same as Reduction A, but 
with 3 arcsec pixels. 



Inverse variance weighted 
flux density summed in 
3 arcsec pixels to produce 
a single flux density map. 
Negative flux density 
from off-chop positions 
also summed in at the 

timestream level(see 

iBorvs et anl2003^ . 



Convolution Form maximum likelihood Same as Reduction A. 


Same as Reduction A, but 


Form maximum likelihood 


point-source fiux density 


with 14.7 arcsec Gaussians. 


point-source flux density 


and noise maps (and S/N 




and noise maps (and hence 


maps) via noise-weighted 




S/N maps) using noise- 


convolution with the differ- 




weighted convolution with 


ential PSF that combines 




a single 14.7 arcsec Gaus- 


14.5 arcsec Gaussians with 




sian (i.e. without the chop 


the chop pattern. 




pattern). 



Cut Pixels No pixel cuts were made 
in fiux density map. Ignore 
pixels with S/N deviating 
more than 4 cr from con- 
volved map. 



No pixel cuts were made 
in flux density map. Ignore 
pixels with a > lOmJy in 
convolved map. 



Remove pixels with < 10 
'hits' in flux density map. 
Ignore pixels with cr > 
5 mjy in convolved map. 



Remove pixels with < 10 
'hits' in fiux density map. 
No pixels were ignored in 
convolved map. 



Source Ex- Positive peaks identified in 
traction the convolved signal map 

above a threshold. A model 
was constructed by cen- 
tering a normalised beam- 
map at the positions of 
the peaks in the con- 
volved map. Normalisation 
coefficients were calculated 
simultaneously, providing 
a minimum noise-weighted 
fit to the unconvolved 

(c) 0000 RAs. M^SB?'^'fi8B: ooo-ooo 



Peaks identified in the 
maximum likelihood S/N 
maps above a 2.5 cr thresh- 
old. Flux densities iden- 
tified at those positions 
in the maximum likelihood 
point source fiux density 
maps. 



Same as Reduction B. 



Same as Reduction B. 
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Figure 1. The 850 fim S/N map with the ofF-beams folded in at the timestream stage (i.e. each chop is treated separately in order to optimise S/N)(left) and corresponding noise map 
(right) of the 485arcmin^ for the LH are shown (specifically for Reduction D). This area only includes Sarcsec pixels that were sampled more than 10 times. On the noise map, 15-arcsec 
diameter circles indicate the positions of the 60 sources in the LH SHADES catalogue. Notice that sources are less often detected in noisier areas of the map; in general they are found 
in areas corresponding to better weather and on a fine scale sources are preferentially found near the corners of the overlapping triangular rows. The SCUBA 8-mJy LH Survey covers 
approximately one quarter of the map in the lower-righthand c orner. Within that r egion is a deep rectangular test strip (p ~ Imjy) that was observed before the rest of the SCUBA 
8-mJy Survey. A similar noise map for Reduction B is found in lMortier et all j200^. 
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Figure 2. The 850 fim S/N map (with the off-beams folded in at the timestream level) (left) and corresponding noise map (right) of the 406 arcmin^ SXDF region are shown (specifically 
for Reduction D). This area only includes Sarcsec pixels that were sampled more than 10 times. On the noise map, 15-arcsec diameter circles indicate the positions of the 60 sources 
in the SXDF SHADES catalogue. In general, sources are more often found on a 'ring' (avoiding a noisier central region taken in predominantly poorer weather), while on a finer scale 
sources are preferentially found near the corners of the overlapping triangular pattern. A similar noise map for Reduction B is found in Morticr ct al. 1,2005.) . 
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Figure 3. Histograms siiowing tlie noise distribution of the 
14.7"FWHM beam-smoothed noise map (pixels with more than 
10 hits from Reduction D) for the LH and SXDF fields in bins of 
A CT = 0.05 mJy beam~^ . The bump in the LH histogram between 
1 and 1.5 mJy arises from the deep strip in the SCUBA 8-mJy 
Survey region. 



3.5 Maps 

The left panels of Figs. 1 and 2 show the point source S/N 
maps for the SHADES fields from Reduction D. The posi- 
tive and negative beams are clearly seen in these maps and 
appear as a result of the differential measurements taken 
at the different chop throws. Noise maps are also presented 
(see right panels of Figs. 1 and 2) with the positions of the 
SHADES catalogue sources marked by circles. Each sam- 
ple of data that is added into a pixel is weighted by the 
variance of the bolometer timestream; the noise map here 
represents the square-root of the variance term. The regular 
grid pattern in the noise maps corresponds to variation in 
the effective observing time as a function of location arising 
from the triangular pointing pattern in the survey and the 
decision to chop the telescope along sky coordinates (RA 
and Dec.) rather than allowing field rotation to smooth the 
observing time over the map. 

Several features in the maps are worth noting. The 
LH map has a non-uniform ch op strategy, since the earlier 
SCUBA 8-mJy Survey region JScott et al.ll2002h was taken 
only with a single 30 arcsec chop (with nod) in Declination. 
This different chop strategy noticeably changes the charac- 
ter of the flux density map; in the SCUBA 8-mJy Survey 
region the noise is clearly spatially correlated in the ver- 
tical direction, and bright sources only have negative side- 
lobes above and below them. In the noise map, however, the 
triangular pattern is repeated across the entire field, since 
the SCUBA 8-mJy Survey was also built from a mosaic of 
overlapping jiggle-maps falling on the same grid adopted for 
SHADES. The deepest region of SHADES is a small rect- 



angular region at the centre of the SCUBA 8-mJy Survey 
where additional data were available from a test strip ob- 
served before the commencement of the full SCUBA 8-mJy 
Survey. Histograms of the noise maps (for Reduction D) are 
given in Fig. |3 The bump in the LH histogram between 1 
and 1.5 mJy arises from the deep strip in the SCUBA 8-mJy 
Survey region. 



4 THE SHADES CATALOGUE 

The SHADES catalogue is compiled using the following 
steps. First, a preliminary joint identification list is com- 
piled by cross-identifying the four independent source lists. 
Consensus SHADES positions and flux densities are then 
determined. Each source flux density is corrected for flux 
boosting (see Section 14.311 . and we reject sources based on 
their deboosted flux density distributions (see Section [4.411 . 
This trimming removes from the preliminary list virtually 
all of the sources which appear to be 2.5-3.5 a in the maps. 
The catalogue we construct here is a robust list, intended 
to be the starting point for follow-up observations leading 
to SED fitting and photometric redshift estimation for indi- 
vidual sources and for the population. 



4.1 Combining partially dependent data 

The four data reductions were carried out independently, 
but since they use the same data, their results are obviously 
not statistically independent. It is relatively straightforward 
to determine if the results of the different reductions are con- 
sistent and we will show that the differences between reduc- 
tions are small compared to the noise in any one reduction. 
Given that the reductions are consistent, it is acceptable to 
choose the result of any single reduction as the final answer, 
but combining the analyses is likely to lead to slightly higher 
precision and reliability. 

The statistics of combining the four reductions is not 
simple, because the degree of statistical independence is 
not well characterised. Differences between reductions would 
arise if systematic errors are present, but might also come 
from random errors and slight differences in weights of the 
input data. Differences in the estimated uncertainties would 
arise if one reduction is genuinely more precise than the oth- 
ers, but also if one reduction over- or under-estimates the 
uncertainties. In combining our reductions we have taken 
a variety of approaches, using the mean value, the median 
or the most precise single reduction, and these choices will 
be described below in the sections on flux densities, source 
positions, and source number densities. We quote as an un- 
certainty the smallest value claimed by a single reduction. In 
doing so we are taking the conservative position that while 
combining the reductions does not decrease uncertainties 
nearly as much as combining fully independent data would, 
it ought not to increase the uncertainty of our most precise 
estimate. 



4.2 Preliminary joint identification list 

An extended source list for each reduction is made of all 
points having S/N > 2.5 in the maps. The S/N threshold 
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4.3 Deboosted flux densities 
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Figure 4. Map flux density minus deboosted flux density versus 
the map-detected flux density for the LH (filled symbols) and 
SXDF (open symbols) SHADES catalogue sources. The line y = 
represents a 1:1 ratio between the flux densities. Map-detected 
source flux densities are deboosted by an average of 1.8 mJy (the 
dashed line), although there are a range of deboosting values, 
because this correction depends on the noise as well as the signal. 



is kept deliberately low to avoid missing genuine sources at 
this stage. 

A preliminary joint list is constructed by identifying 
sources for the four extended lists which are within 10 arcsec 
and which are seen with S/N > 3 in at least two reductions. 
This preliminary list contains 94 source candidates in the 
SXDF and 87 in LH. 

We compute a SHADES map flux density for each 
source using the following recipe. We compute a raw flux 
density likelihood distribution by multiplying the individ- 
ual Gaussians constructed from the flux and noise in each 
reduction and normalise the resulting distribution. We take 
the SHADES map flux density to be the maximum likeli- 
hood. We use the lowest quoted error as the 1 a SHADES 
map noise uncertainty. Because the data are common to the 
four reductions, adding errors as inverse variances would se- 
riously underestimate the net uncertainty, while simply aver- 
aging the uncertainties allows an imprecise single estimate to 
lower the combined uncertainty, which does not make sense. 
Reductions B and D agree well on flux densities in both 
flelds and also claim the smallest photometric uncertainty, 
so data from those two reductions dominate the weighted 
mean flux density for those sources which they extract. We 
find that the deviations between any single reduction and 
the group flux density are smaller than the adopted mea- 
surement noise, as expected. The details are in Section [4.61 
The ratio of SHADES map flux to minimum noise is listed 
as S/N in Table H 



We have empl oyed 
ICoppin et alj Jic 



the simple Bayesian recipe of 
to correct the preliminary source list 



for effects of flux density boosting in submillimetre maps. 
Details of how we deboost flux densities from the joint can- 
didate list are given below. After deboosting, we are left with 
a catalogue of 120 robust SHADES sources, coincidentally 
60 sources per field. 

The submillimetre source count density in the flux range 
of the SHADES survey falls very rapidly compared to the 
width of the approximately Gaussian noise distribution of 
our maps. Therefore we expect to identify an excess of 
low flux density sources whose locations happen to coin- 
cide with positive noise and whose apparent flux densities 
have therefore been increased above the survey's S/N limit. 
This 'flux boosting' is a well-known effect in low S/N sub- 
millimetre maps and must be quantified so that individ- 
ual source flux densities can be corrected to represent the 
best estimate of the true underlying flux density of each 
source. This boosting is sometimes referred to as Malmquist 
bias, although more properly that refers to seeing sources 
more luminous than average in a magnitude-limited sur- 
vey because of a spread in luminosity, i.e. there is still 
Malmquist bias even when there is no measurement error 
llHendrv. Simmons fc Newsamlll993t) . Flux density boosting 
is also distinct from Eddington bias, which properly refers to 
the effect on the number counts rather than the individual 
source flux densities. 

Each candidate sour c e's flu x density is deboosted fol- 
lowing the ICoppin et alJ (|200^ Bayesian recipe, resulting 
in a posterior flux density probability distribution which is 
altered from the Gaussian probability distribution inferred 
from the maps. The probability distribution for an individ- 
ual source's intrinsic flux density Si given that it is found at 
the observed flux density So ± CTo, is p{Si\So, cto), which is 
factored using Bayes' theorem: 



p(gi)p(5'o,ao|5i) 



(1) 



This expression states that the intrinsic flux density distri- 
bution of the source is the likelihood of observing the data 
(right hand term in the numerator) weighted by the prior 
distribution of flux densities (left hand term in the numer- 
ator). A Gaussian noise distribution is used for the likeli- 
hood of the data (the photometric error distribution). We 
use an informative prior for pixel flux densities constructed 
from existing knowledge of the 850 /im extra-galactic source 
counts and the SHADES observing strategy. This inherently 
gives more weight to lower flux density sources, since they 
are more numerous. To co nstruct the p r ior, w e flt the broken 
power- law model given in lScott et al.l ll2002l) to the counts 
of Borvs ct al. (2003) and constrained by the lensing counts 
at fainter fluxes, but have checked that mild deviations from 
this make no significant difference. The model is of the form 



dN 

Us 



S" 



s_Y fs_ 

S' ^ I S" 



(2) 



Artificial skies are generated consistent with Equation |5| 
Poisson statistics, random source locations and each sky 
is sampled with the SHADES chopping pattern and point 
spread function to obtain a pixel flux density distribution. 
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Because our maps are differential with zero mean, and 
because tfiey contain noise, tfie expression in Equation^can 
have a negative tail. The fraction of the posterior distribu- 
tion p(Si\So,<7o) having 5*1 < is taken as the probability 
that a given source is falsely detected. 

The idea of using Bayes' theorem to find a posterior 
estimate of the flux density using the source counts as 
a prior app e ars to have been first clearly written down 
by IJauncevI lll968l) . as a way to correct survey-detected 
ra dio sources for flu x density biasing (as first pointed out 
bv lEddingtonllml. Other papers which discus s similar 
ideas^ include llVIurdoch. Crawford fc Jauncev i WTZ ), 
Schmidt fc Maccacard Jl986D . iHogg fc Turneil Jl99d) . 



Wang. Cowie fc Barged t004) . and iTeerikorpil J200. 



On average, the deboosting reduces the source flux den- 
sity by ~ 1.8 mjy (see Fig. 01 , increases the width of the pho- 
tometric error distribution by about 10 per cent, and ren- 
ders the shape of the resulting distribution to be skewed and 
non-Gaussian. The details of these effects depend both on 
the observed signal. So, and the observed noise, CTo, and not 
just on the S/N. The effects are larger for sources extracted 
from noisier regions of the maps. We have plotted two ex- 
amples of posterior distributions taken from the LH region 
in Fig. 13 The first is a bright source and the other is dim; 
one readily sees that the skew of the distribution is more 
pronounced in the dim case. We note that pr e vious sub- 
millimetre survey s such as those of IScott et al.l ll2002f) and 
lEales et al.l (l2000t) assessed the level of flux density boosting 
through simulations by comparing input source flux densi- 
ties to extracted source flux densities for a 3 a catalogue 
cut and noi se map RM S valu es of ~ 2.5 mJy and 1 mJy, 
respectively. IScott et alJ (|20o3) found a boosting factor of 
about 15 per cent at SmJy and 10 per cent at 11 mJy, while 
lEales et al.l i2000ri quote a median boost factor of 1.44 with 
a large scatter: these correction factors were detemined as 
a function of flux density only. In SHADES, each source is 
individually corrected for flux density boosting by a factor 
determined as a function of the map-detected flux density 
as well as the noise estimate. 

The deboosting recipe has been successfully tested 
against follow-up p hotometry for individual sources in 
ICoppin et al.l (|2003). Its performance in returning the in- 
put source distribution is tested in Section fS. 1.1 1 

Since some of the LH data were taken originally for 
the SCUBA 8-mJy Survey llScott et alJl2002^ using a sin- 
gle chop throw, the effect of using this observing scheme on 
the reported deboosted flux densities was investigated by re- 
calculating the prior and flux density posterior probability 
distributions. Reported flux densities were different for less 
than 10 per cent of the SHADES sources, and at most by a 
negligible ±0.1 mJy. In addition, 3 extra sources near the re- 
jection threshold made it into the flnal LH source catalogue. 
For simplicity, we used the prior based on the multiple chop 
SHADES observing scheme to deboost all of the sources, 
since it is representative of the majority of the data. 

The possible effects of clustering have not been in- 
cluded in creating the prior distribution used in deboost- 
ing fluxes. We have checked that clustering at the levels 



In stellar astronomy the sa me effe ct on parallax measurements 
is known as Lutz-Kelker biaslLutz fc Kclkcr L19I2.). 



anticipated for SHADES sources has a negligible effect on 
the deboosted flux density distributions. Using 50 realisa- 
ti ons of the phenomenologic al galaxy formation model used 
in Ivan Kampe n et alj (|200^, with a clustering strength of 
6o — 10 arcsec, we create a noiseless distribution of map 
pixel fluxes and use it to construct a new prior distribution 
and deboost sources with similar flux densities and S/N as 
the SHADES sources. We compare the posterior flux den- 
sity distributions to those calculated using a prior with the 
same input source count model, but with randomised po- 
sitions, i.e. no clustering. We find negligible differences be- 
tween the distributions. A larger effect on the shape of the 
posterior fiux density probability distribution comes from 
using a much steeper source count model at the faint end of 
the number counts, which has a small but noticeable effect 
on the shape of the posterior fiux density probability distri- 
bution at flux densities below 1 mJy, making sources more 
likely to pass the catalogue cut. We therefore claim to be 
making the most conservative catalogue cut, as we use less 
steep source counts at the faint end that we believe are a 
balance between faint lensing and bright blank-fleld counts. 



4.4 850 /im catalogue membership 

The flnal cut in catalogue membership is the requirement 
that each accepted source has less than 5 per cent of its 
posterior probability distribution below OmJy, or P{Si < 
OmJy) < 5 per cent. This threshold is a good balance be- 
tween detecting sources while keeping the number of spuri- 
ous detections to a minimum. 

One could tune the catalogue membership using the 
thresholdin g technique called the False Dis covery Rate 
(FDR ; see iBeniamini fc Hochber3 Il995l and iMiller et all 
120011) in order to control the average fraction of spurious 
sources in the catalogue. In Fig.|S]we have plotted the pos- 
terior null probability for each source candidate (i.e. the per- 
centage of the posterior flux density probability distribution 
which is below OmJy) ranked in ascending order. This plot 
illustrates how our choice of probability cut-off for individ- 
ual sources occurs comfortably before the regime where the 
FDR increases dramatically. It is also coincidentally where 
the number of sources times the FDR approximately equals 
1. These plots show that as one pushes beyond ~ 60, the 
number of spurious sources becomes ^ 1. 

Given the number of beams in each map, one expects 
approximately flve 3 a peaks at random. If we relaxed the 
null probability cut we would increase the number of sources, 
but the chances of random noise peaks making it into the 
catalogue would then become important. The average of the 
null probabilities is 1.5 (2.0) per cent in the LH (SXDF), and 
this can be interpreted as an overall FDR for the catalogue. 

Effective flux density and noise cuts in each catalogue 
are shown in Fig.|7] No sources with observed S/N < 3.2 are 
kept in the flnal catalogues, with the majority of the detec- 
tions lying above 3.5 a (see Fig. |H|for the S/N distribution 
of the SHADES catalogue sources). 

The flux densities we quote in Table |21 are median flux 
density estimates and the quoted errors correspond to the 
central 68 per cent of the posterior flux density distribution. 

We remark that there are substantially more high S/N 
sources with null probabilities below 5 per cent (i.e accepted 
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Flux Density [mjy] Flux Density [niJy] 

Figure 5. Posterior flux density probability density functions calculated from Reduction D by the method described in Section 14.31 are 
shown as dot-dashed curves for a high S/N source, LOCK850.5, and a low S/N one, LOCK850.52. The thick dashed vertical line indicates 
the average SHADES map flux density before deboosting. Notice how asymmetric the low S/N posterior flux density distribution is. 
For comparison, distributions calculated by alternate methods from Reductions B and C are shown for the same sources. Each group's 
map flux density is within 0.5 mJy of this average flux density before deboosting. In these cases, the SHADES posterior flux density 
distribution (not shown) follows closely the shape of Reduction D's posterior flux density distribution. All of the distributions shown 
have been normalised to have unit area. The distributions from Reduction C are lower resolution due to the choice of bin size for the 
simulations used to determine the posterior distribution and are truncated below 2 m Jy (see Section 15.21 . Posterior distributions from 
Reduction B are truncated outside the region between 0.5-2 times the map flux density. 




Observed Flux Density |mJy| Observed Flux Density ImJyl 



Figure 7. Effective cuts in the LH (left) and SXDF (right) source catalogues are shown. For all candidate sources, the SHADES observed 
S/N is plotted against observed flux density as circles (the size of a circle is proportional to the map noise level where a source was found). 
The dashed lines show observed noise levels of 2, 3, and 4 mJy. Sources are retained in the catalogue if the total posterior probability that 
the flux density is zero is less than 5 per cent (filled circles) . Only the sources lying above the solid curve satisfy this criterion. The effect 
of the rapidly falling source count spectrum is visible as the curvature of this cut-off, i.e. the rise in the required S/N with increasing 
noise. An apparent 3.5 cr source in a noisy region of the map has a higher flux density than a 3.5 cr source found in a quiet region; bright 
sources are rare, so this source is less likely to be genuine than a quieter dim source is. Notice that there are dramatically fewer sources 
detected at flux densities above 10 mJy in the SXDF as compared to the LH field, even though the number of fainter sources is similar. 
Also note that essentially none of the SXDF detected sources have noise in the 1— 1.5mjy range, as compared to the LH, which includes 
the lower noise SCUBA 8-mJy Survey region. 



© 0000 RAS, MNRAS 000, 000-000 



12 Coppin et al. 




Number of Sources 850)im map S/N for the SHADES Sources 



Figure 6. Percentage likelihood of zero flux density. This plot 
shows the posterior probability that each source has a flux den- 
sity < Omjy. This probability is plotted for each source candidate 
in the LH (filled symbols) and SXDF (open symbols) in ascending 
order. Only those candidate sources for which the null probability 
is less than 5 per cent (corresponding to locations in the figure 
below the horizontal dashed line), are kept in the SHADES cat- 
alogue. In each of the LH and SXDF there are 60 such sources. 
Notice the comparatively small number of sources in the SXDF 
with very low null probabilities compared with the LH. Notice 
also that if we had chosen a cut below 5 per cent there would 
have been more LH than SXDF sources, while the opposite would 
have been true if we had chosen a cut above 5 per cent. Given 
the number of beams in each map, one expects typically five 3 cr 
peaks at random, and perhaps three of these survive the deboost- 
ing process. The solid line shows the percentage of the catalogue 
comprised of such sources as a function of the number of sources 
in the catalogue, P = 3/N, which coincidentally crosses in the 
same place. 



Figure 8. Distribution of the map-detected S/N ratios (SHADES 
catalogue values) for the 60 LH and 60 SXDF sources in bins of 
A S/N = 0.5. We note that all of the SHADES catalogue sources 
were initially detected above a S/N of 3.2, whereas the prelimi- 
nary list contained sources > 2.5. 



in the SHADES catalogue) in the SHADES maps of the LH 
region than in the SXDF (see Fig. IHJ. 



4.5 Final SHADES catalogue 

The SHADES catalogue is given in Table |21 Gaps in the 
source numbering sequence indicate sources that were re- 
jected from the preliminary catalogue because they either 
failed to be detected by at least two groups with S/N > 3, 
or because P{Si < OmJy) > 5 per cent. Comments on par- 
ticular sources are noted in Appendix 

For the first time in a submillimetre-selected survey, a 
careful estimate of the unbiased flux density of each source 
is provided (cf. Section [4.31 . 

Table |21 also contains 3 a upper limits for the 450 /im 
flux densities of each SHADES source (see Appendix IB3. 31 
for details). 
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Name (lAU) 


Nickname 


RA 


Dec. 




^850 


map S/N 


map 




5450 


Other IDs 








(J2000) 


(J2000) 




(mjy) 




(m 


Jy) 


(mJy) 


/Notes 


SHADES 


J105201+572443 


LOCK850.1 


10''52"' 


01? 42 


57° 24' 43' 





8.85 {±;:;j) 


8.6 


1 


1 


< 47« 


LE850.1,LE1100.14,t;?Xs 


SHADES 


J105257+572105 


LOCK8B0.2 


10''52"' 


57! 32 


57°21'05' 


8 


13.45 (±1;}) 


6.8 


2 


1 


< 123 


LB1100.1,<7N 


SHADES 


J105238+572436 


LOCK850.3 


10^52™ 


38! 25 


57°24'36^ 


5 


10.95 (±}:|) 


6.4 


1 


9 


< 34 


LB850.2,LE1100.8,(? 


SHADES 


J105204+572658 


LOCK850.4 


10''52"' 


04! 17 


57°26'5S^ 


9 


10.65 (±;:|) 


6.4 


1 


8 


< 134 


LE860.14,<:5 


SHADES 


,7105302 + 571827 


LOCK850.5 


10''53"' 


02! 62 


57° 18' 27' 





8.15 (±2°) 


4.9 


2 





< 107 




SHADES 


,7105204+572526 


LOCK850.(i 


lo''52"' 


04! 13 


57° 25' 26^ 


3 


6.85 (±i;|) 


5.8 


1 


3 


< 77 


LE860.4,<? 


SHADES 


,7 1 ( J 5 3 (J 1 + 5 7 2 5 5 4 


LOCK850.7 


lt)^^53™ 


0l!40 


57° 25' 54' 


2 


8.55 (±}|) 


5.3 


1 


8 


< 85 




SHADES 


,7105153 + 571839 


LOCK850.8 


lo''51°' 


53! 86 


57° 18'39' 


8 


5.45 (±11) 


5.2 


1 


2 


< 31 


LE850.27,<:?Xt 


SHADES 


J105216+572504 


LOCK850.9 


10''52"' 


16! 09 


57°25'04^ 


1 


5.95 (±l\e) 


4.7 


1 


6 


< 68 


LE860.29,^ 


SHADES 


J 105248+573258 


LOCK8B0.10 


10^52™ 


48? 61 


57° 32' 68" 


6 


9.15 (±11) 


4.6 


2 


6 


< 366 


t 


SHADES 


J105129+572405 


LOCK850.il 


10^51™ 


29? 53 


57° 24' 06" 


2 


6.25 (±1:1) 


4.6 


1 


7 


< 62 




SHADES 


J105227+572513 


LOCK850.12 


10^52"^ 


27?61 


57° 25' 13" 


1 


6.15 (±1:1) 


4.6 


1 


6 


< 35 


LE850.16,LE1100.16,^Xt 


SHADES 


J105132+573134 


LOCK850.13 


10''51'" 


32? 33 


57°3l'34V 


8 


5.65 (±|;|) 


3.9 


2 





< 109 


t* 


SHADES 


J105230+572215 


LOCK850.14 


10''52'" 


30? 11 


57°22'15( 


6 


7.25 {±J;g) 


4.8 


1 


8 


< 96 


LE850.6,LE1100.5,^ 


SHADES 


J105319+572110 


LOCK850.15 


10''53"' 


19? 20 


57° 21' 10' 


6 


13.25 (±i:l) 


4.6 


3 


8 


< 149 


^Xt 


SHADES 


,7105151+572637 


LOCK850.16 


lo'^si" 


5l!45 


57°26'37! 





5.85 (±1:1) 


4.3 


1 


7 


< 67o 


LE850.7,<:? 


SHADES 


J105158 + 571800 


LOCK850.17 


10''51"' 


58! 25 


57°18'00^ 


8 


4.75 (±1:1) 


4.6 


1 


3 


< 55 


LE850.3,^ 


SHADES 


.7105227+572217 


LOCK850.18 


lo''52"' 


27! 69 


57°22'l7^ 


8 


6.05 {±^;?) 


4.3 


1 


8 


< 84 


SDS.20,t 


SHADES 


.7105235 + 5 73119 


LOCK850.19 


10^^52™ 


35! 71 


57° 31' 19' 


1 


5.15 (±2°) 


3.9 


1 


8 


< 87 




SHADES 


.1105256 + 573038 


LC)CK850.21 


10''52°' 


56! 86 


57°3t)'38' 


1 


4.15 (±2'°) 


3.6 


1 


7 


< 70 




SHADES 


.1105137+573323 


LOCK850.22 


10''51"' 


37! 55 


57° 33' 23' 


3 


7.00 {±i:i) 


4.0 


2 


7 


< 76 




SHADES 


J105213+573154 


LOCK850.23 


10''52"' 


13! 74 


57°3l'54( 


1 


4.35 (±t'l) 


3.7 


1 


7 


< 60 




SHADES 


J105200+572038 


LOCK850.24 


10''52'" 


00? 23 


57° 20' 38" 


1 


2.75 (±il) 


3.6 


1 


1 


< 32 


LE850.32,<? 


SHADES 


J105240+572312 


LOCK850.26 


10** 52™ 


40? 9 5 


57° 23' 12" 





5.85 (±1;*) 


3.9 


2 


1 


< 48 


t 


SHADES 


J105203+571813 


LOCK850.27 


10''52'" 


03?57 


57° 18' 13" 


6 


5.05 (±1:1) 
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Tabic 2: The 850 fim SHADES catalogue for the LH and SXDF regions. Uncertainties 
in physical positions are 3.2 arcsec in RA and 3.2 arcsec in Dec. {the same uncertainties 
have been adopted for all SHADES sources in both fields). Estimates of the true unbiased 
median flux density of each source is given, with accompanying error bars representing the 
68 per cent confidence bounds of the (non-Gaussian) deboosted flux density distribution 
(cf. Section for the Bayesian estimate of deboosted flux density). We also provide 

the combined map S/N and noise estimates (i.e. values prior to deboosting). We do not 
claim any 450 /xm detections of SHADES sources; 3 (T upper limits limits are given i n the 
penultimate column. Correspon ding SCUBA 8-mJy Survey IDs ISc^^^^^!h| |20o3) and 
new IDs from the reduction of ^^^^^^^^^^^^^^^^^^^^^ ^^^^ are listed in the final 
column for reference, along with an indication of how many groups detected each source. 
See Appendix^^for detailed notes on some of the sources. 

Legend: 

's? This source was identified by all four groups with a S/N > 3. 

I This source was identified by three groups with a S/N > 3. 

X This source was identified by two groups with a S/N > 3. 

* This source was identified by one additional group with 2.5 < S/N < 3. 

** This source was identified by two additional groups with 2.5 < S/N < 3. 

K This source is mildly affected by the noise spike (cf. Appendix^^. 

Xsi Xt This source has a relatively poor spatial or temporal x 

O This source has a hint of flux at 450 fim measured at the SHADES catalogue positions 
at a significance level > 3 cr (see Fig.^^^. 
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^.5.1 Comparison with the SCUBA 8-mJy Survey 

The new SHADES catalogue was cross-match ed with the 
SCU BA 8-mJy Survey > 3.0 (t source catalogue iScott et alJ 
|20o3), which used a subset of our data. We failed to re- 
detect 2/12 of the >4cr sources (LE850.5 and LE850.9), 
5/9 sources with published S/N between 3.5 and 4.0 
(LE850.13, LE850.15, LE850.17, LE850.19, and LE850.20), 
and (not surprisingly) 12/15 sources with S/N between 3.0 
and 3.5 (LE850. 22-26, LE850.28, LE850.30, LE850.31, and 
LE850. 33-36 ) . The s e findings ar e similar to t hose given in 
llvisoii et ai] J2002h. iMort ier et'al. ( J0Q5) and llvison et all 
l|2005^ . llvison et a l. ( 2002) rejected "lE850.9, LE850.10, 
LE850.15, and LE850.20 due to the lack of associated ra- 
dio counterparts, combined with the fact that they were 
found in noisy regions of the map (ct > 3 mjy). These sources 
are also rejected in our analysis, except in the case of 
LE850.10, since this source is re-detected in the SHADES 
data (LOCK850.60), albeit with a lower S/N than that 
found in the SCUBA 8-mJy Survey (~ 3.4 ct as compared 
with the previous 4.2 ct detection). We find that there is less 
than a 4 per cent chance of LE850.10 having a true flux 
density of < mJy and therefore it survives the deboosting 
cut. This source also has a tentative radio identification (Ivi- 
son et al. in preparation). See Table|21for the corresponding 
new SHADES measurements of the SCUBA 8-mJy Survey 
sources. 

The SHADES catalogue was also cross-matched with 
the re-red uction of the SCUBA 8-mJy Survey > 3.0(7 source 
catalogue iScott. D unlop & Scrieant 2006), which includes 
some additional data and i mprovements made to the reduc- 
tion methods. In summary, IScott. Dunlop fc Serjeant (1200611 
re-detected all of the 36 original SCUBA 8-mJy Survey 
> 3.0 o sources (though 4 sources originally detected above 
3.0 a dropped down to below 3.0 a; LE850.25, LE850.29, 
LE850.30, and LE850.31), and found 8 new 3.0 < S/N < 
3.7 sources. SHADES re-detects sources SDS.16 (3.7 a) 
and SDS.20 (3.4 a), but fails to re-detect SDS.25, SDS.32, 
SDS.33, SDS.36, SDS.38, and SDS.40 (new < 3.4 a sources). 
See Table 121 for the corres ponding SHADES measurements 
of the sources found in the lScott. Dunlop fc Serieanti (l2006l) 
re-reduction of the SCUBA 8-mJy Survey data. 



4.6 Flux density comparison 

We checked for systematic effects between the reduction 
group flux densities in order to quantify their effect on the 
adopted SHADES map flux density. We limit the discus- 
sion that follows to a subset of sources that all groups 
flnd at > 2.5 a: 54 in the LH; and 58 in the SXDF. Un- 
like in the astrometry comparison (see later Section 14.711 . 
we have no 'true' flux density to compare the mean flux 
densities against, so only an inter-group comparison can be 
performed. 

Flux density comparison scatter plots were produced, 
such as Fig. ^ for each pair of groups. We noticed immedi- 
ately that one reduction showed a systematically lower flux 
density than the other groups that appeared to depend on 
the pointing strategy of the observations. After finding and 
correcting this error, the scatter in source fiux densities be- 
tween all groups appeared small on average, with no system- 
atic offset apparent in any one group compared to another 



(except at the high flux density end, where the photometric 
errors are also very large). We feel it is important to mention 
explicitly those very few occasions when the results of our 
comparisons are used to guide the details of any reduction, 
so that a critical reader will understand how close to a blind 
analysis the SHADES reductions are. 

One might expect the choice of sky opacity correction 
factors or FCFs (i.e. using monthly FCFs versus nightly 
measurements) to come into play at about the 2 per cent 
level for such low S/N data taken in more or less uniform 
weather conditions. The RMS scatter between flux densities 
reported for LH sources by Reductions B and D is 0.9 mJy 
(see Fig. ^ and is similar between the other groups. Similar 
results were found for the SXDF. We find that the photom- 
etry errors we might have introduced through differences in 
judgment are 1/3 as large as the total uncertainty in flux 
density, which is noise-based. We therefore claim that the 
gain differences are small (i.e. less than 5 per cent) and 
therefore unimportant for these low S/N data (they become 
important for ~ 10 a sources, of which we do not find any 
in this survey). This demonstrates that we are extracting 
flux densities well. Using monthly averaged FCFs (Reduc- 
tion D) versus using nightly measurements (Reductions A, 
B and C) appears to give indistinguishable answers, and 
therefore the systematic error introduced from the former 
technique and the instantaneous measurement uncertainty 
in the latter technique are both insignificant with respect to 
the photometric errors in our survey. 

4.7 Astrometry comparison 

Intercomparing the positions found by the four independent 
reductions allows us to check for systematic errors and es- 
timate the uncertainty of our astrometry. When comparing 
flux densities, differences between reductions can be mea- 
sured but the 'true' 850 ^m flux density is not known, so 
systematic errors may be difficult to isolate. The situation 
is much simpler with positions. When sources have clearly 
identified radio counterparts the precision of the positions 
determined from the radio data is much higher than that 
available from the 850 yLxn data. Effectively, one knows the 
'true underlying position' for these sources and sensitive 
tests for systematic errors are possible. In this section we 
analyse the positions of a subset of the sources which have 
clear and compact radio counterparts determined by Ivi- 
son et al. (in preparation) and which have been identified 
at > 2.5 a in all four reductions. These criteria, designed to 
facilitate clear comparison of the four reductions, yield 17 
and 24 sources in the LH and SXDF, respectively. The full 
analysis of the alignment of SHADES sources with the corre- 
sponding radio data is made in Ivison et al. (in preparation) . 

Upon initial comparison it was clear that the positions 
of sources in the LH determined by one reduction were dis- 
placed to positive RA by just over 2 arcsec compared to posi- 
tions determined by any of the other reductions or compared 
to the radio positions. This systematic effect was traced to 
errors in the use of the hastrom routine in idl-astrolib and 
has been corrected. Other than the correction of this small 
error, nothing has been adjusted to bring the reductions into 
agreement with each other or with the radio-determined po- 
sitions. 

The positions determined by one reduction (arbitrar- 
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Figure 9. The diflfcrence of the 850 fira flux density as measured by Reductions B and D is plotted against the flux density measured by 
Reduction B for the LH sources that were found in all four reductions (left panel), and a histogram of flux density ratios of Reductions 
B and D (i.e. B/D) for the LH (right panel). In the left panel the horizontal solid line indicates a ratio of unity, while the dashed lines 
show the 1 a scatter expected if each of the two reductions had independent photometry errors corresponding to a S/N of 3.5. The scatter 
between these two reductions is substantially smaller than photometry errors, indicating that systematic errors associated with data 
reduction choices do not substantially influence our final flux density measurement values and uncertainties. Reductions B and D are 
selected for comparison here because their calibration procedures differ most among the four reductions; any other pair of reductions is 
likely to demonstrate at least this level of agreement. In the right panel the smooth (dashed) curve is the flux density ratio distribution 
expected for two independent measurements with the approximate limiting SHADES S/N (3.5(t). The narrowness of the histogram 
compared to the relative width of the Gaussian demonstrates that systematic errors in calibration are small compared to photometric 
uncertainties for all sources in the SHADES catalogue. The two sets of flux density values are consistent with each other, while the best 
fit mean flux density difference is 4 per cent. Since the uncertainty of the mean of all 60 sources in the SHADES LH catalogue is also 4 
per cent, systematic fiux density errors associated with differences in data reduction strategy are completely unimportant for individual 
sources. 



ily, Reduction D) in comparison to the mean 850 /im po- 
sition for all the sources in the SXDF are plotted in the 
left hand panel of Fig. [TOl The sources with compact radio 
IDs from Ivison et al. (in preparation) are shown as plusses, 
while other sources in the SXDF are shown as diamonds. 
Both sub-samples have similar means and distributions, so 
we conclude that the analysis of positions for the restricted 
sub-sample provides a description of astrometry errors which 
is valid for the full list. Table^lists the RMS displacements 
of source positions determined by each reduction from the 
unweighted mean position determined by the four reduc- 
tions, separately for the two SHADES fields. Pixelisation 
inevitably adds d/^/vl where d is the pixel size, in quadra- 
ture to the RMS astrometry errors^. This is 0.9arcsec for 
3 arcsec pixels and has not been subtracted from the data in 
Table|3 Even so. Reductions C and D, with 3 arcsec pixels, 
have RMS displacements from the 850 /im mean which are, 
if anything, lower than the displacements of the reductions 
using smaller pixels. It is perhaps not surprising that using 
small pixels (smaller than ~ ^^^hm least) in reconstruc- 
tion of the data does not appear to add any astrometric 
precision (cf. ICflnH^BEll)- 

The right hand panel of Fig. 1101 shows the devia- 
tion of the mean offset of the position determined by the 
four SHADES reductions relative to the presumably correct 

* The location of a detected source is uniform between —1/2 
and 1/2 of the centre of a pixel of size 1. The variance of the 
error made by quoting the pixel centre is the mean-square error 
= S%x-^dx = l/l2. 



radio-determined positions from Ivison et al. (in prepara- 
tion) . Notice that the scatter is much larger than the scatter 
in the left hand panel. We conclude that the four reductions 
are accurately extracting the location of the peak in the sub- 
millimetre flux density from the data, and that this peak is 
displaced from the true source location due to noise in the 
850 ^m data. In a careful comparison, Ivison et al. (in prepa- 
ration) confirm that the offsets scale as expected with the 
submillimetre beam size and S/N. We quote a positional 
uncertainty in RA and Dec. offsets of 3.2 arcsec and find no 
evidence for an overall mean astrometric error. It is inter- 
esting to note that the unweighted mean position from the 
four reductions is a better predictor of radio position than 
is obtained from any single reduction. Taken together the 
panels in Figure [111 indicate that the reductions do a good 
job of determining the positions implied by the submillime- 
tre data, but that those positions have a few arcsec scatter 
with respect to the true underlying source positions. 



4.8 Source deblending 

The source extraction methods of Reduction C (D) are in- 
sensitive to finding sources closer than 10 (18) arcsec. There 
is thus a potential problem with source blending, since the 
detection of very near neighbours would be evidence for clus- 
tering of SMGs. 

Motivated by the discovery of two additional sources 
near two bright SMGs dete cted in the GOODS-N SCUBA 
map by IPope et alJ (l2005l) . double Gaussians with vari- 
able positions and amplitudes were fitted to the 5 bright- 



© 0000 RAS, MNRAS 000, 000-000 



SHADES suhmillimetre maps, catalogue and number counts 17 



Table 3. Astrometry Precision: The submillimetre positions of a subset of SHADES sources are compared to each other and to the 
radio positions found in Ivison et al. (in preparation). The columns marked 'Std. dev. w.r.t. 850 /jm' show the standard deviations in 
arcseconds of each reduction with respect to the unweighted mean of all four reductions. The column marked 'Std. dev. w.r.t. Radio' 
shows the standard deviations in arcseconds with respect to positions determined by Ivison et al. (in preparation). Columns marked 'Net' 
list the quadrature sums of the RA and Dec values across both fields. The positions of the radio sources are known to a higher accuracy 
than is possible using the SCUBA data alone, assuming that all radio counterpart IDs are correct. The pattern of variances in the table 
indicates that the peak in the submillimetre data is displaced from the true source position due to noise, but that all reductions find 
consistently the same displaced peak location because the submillimetre noise is common to all four reductions. The unweighted mean 
position of the four reductions is a better predictor of radio position than any single reduction. 



Std. dev. w.r.t. 850 fiia Std. dev. w.r.t. Radio 



SHADES Lockman SXDF Net 850 Net Radio 

Reduction RA Dec. RA Dec. RMS RMS 



A 2.33 1.94 1.80 2.91 2.29 3.46 

B 1.13 1.20 1.40 2.23 1.55 3.06 

C 1.26 1.14 1.49 1.71 1.42 2.77 

D 0.93 1.43 1.38 1.36 1.29 3.01 
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Figure 10. The left panel shows the offset of the source positions inferred by Reduction D from the mean of the positions found 
by all four reductions for all sources in the SXDF which are found in every reduction. The plusses show those points which are also 
included in the comparison with radio positions (see right panel and Table 1^ . The diamonds are the remaining sources and they do 
not appear to be distributed very differently than those we have compared to the radio counterparts. The dashed ellipse indicates \/2 
times the one-dimensional variances in RA and Dec; it contains ~ 2/3 of the sources, as expected. The right panel shows the offset in 
position of the mean of the four SHADES reductions relative to the radio data of Ivison et al. (in preparation) for a subset of sources 
which have clear well defined radio positions. The radio data are collected at substantially higher angular resolution and S/N, so the 
scatter here is presumed to be dominated by scatter in the submillimetre data. As in the left panel, the semi-major and semi-minor 
axes of the dashed ellipse are \/2 times the positional variances in RA and Dec. This error ellipse should contain ~ 2/3 of the sources. 
The vertical and horizontal dashed lines show the mean displacements of SHADES locations from the radio; these differences are not 
statistically significant. Taken together these panels indicate that the reductions do a good job of determining the positions implied by 
the submillimetre data, but that those positions have a few arcsec scatter with respect to the true underlying source positions. See Ivison 
et al. (in preparation) for a more detailed comparison. 
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est sources in each field found by Reduction D. However, 
no convincing evidence was found for favouring two sources 
over one. The '^^^ typically lower when fitting double 
Gaussians, but the second source was never bright enough 
to be classified as a detection under our criteria. We note 
here (and in Appendix (CJ that SXDF850.1 appears quite 
extended in the NS direction, although the is not lower 
for two sources than for one. 



5 RESULTS: DIFFERENTIAL SOURCE 
COUNTS 

An important quantity which can be derived from the data 
is an estimate of the number of sources as a function of fiux 
density. SHADES is a single, uniform survey which has ap- 
proximately doubled the total area of all SCUBA-observed 
blank-field surveys. While the dynamic range in fiux den- 
sities is not as broad as is available from a compilation of 
data which includes the deepest surveys (particularly those 
associated with foreground gravitational lenses), the results 
here are the most robust obtained so far in the range of 
fiux densities where SHADES is sensitive. In this section we 
provide reliable estimates of the differential 850 /im number 
counts for the first time. Difi'erential counts ofi'er an advan- 
tage compared to integral source counts, since each estimate 
of the number of sources in a flux density bin does not de- 
pend on the counts at brighter fiux densities and thus will be 
much less correlated. This makes fitting the results to models 
of source counts more straightforward. Nevertheless, in the 
following section, |S| we integrate the differential counts to 
estimate the cumulative source distribution for comparison 
to previous data. 

Although all of the sources in the SHADES catalogue 
are at least 3 cr in the maps, and therefore have raw measured 
fluxes which are typically above 6 mjy, the well known steep 
submillimetre source flux density distribution implies that a 
source which is detected at 6-8 mJy is as likely to be a 3- 
4 mJy source accidentally observed with a positive noise fluc- 
tuation as it is to be a genuine 7mJy source (i.e. Malmquist- 
like bias; cf. Section 14.311 . The number of sources we have 
observed depends on the number density of sources down 
to quite faint flux densities, well below our nominal 6mJy 
limit. Fits of differential source counts to our data will there- 
fore provide constraints on the number of sources per square 
degree starting at about 3mJy. 

To obtain estimates of differential source counts one 
must estimate completeness, flux density bias, survey area, 
spurious detection rates, and Poisson counting errors. We 
have developed two analysis approaches to perform these 
tasks. One of these is a 'direct estimate', which works with a 
list of sources and their deboosted fluxes and sums the asso- 
ciated probability densities to obtain the parent source den- 
sity spectrum. The second method, 'parametric fitting', self- 
consistently estimates the prior source density spectrum, the 
FDR and the source deboosting. Direct counts use a fixed 
informative prior to perform the deboosting. This is a max- 
imal use of existing information from past SCUBA surveys. 
The parametric approach leaves the prior as a free param- 
eter. This second technique is probably more conservative, 
and the amount of noise in the answer strongly depends on 
how the prior is parameterised. 



Although the SHADES catalogue is extremely robust 
in terms of a low expected false detection rate, it has a 
complicated selection function and is not necessarily opti- 
mal for measuring the source counts. Because it is based 
on four different analyses the selection criteria are hard to 
Monte Carlo, so we do not use the SHADES catalogue to 
derive constraints on 850 /im source counts. In particular, 
one would like to use statistical information from sources 
which may not individually be detected with high signifi- 
cance. Therefore, source count spectra have been determined 
independently from the provisional source lists arising from 
each reduction. Variations of the direct estimation approach 
have been applied to data from Reductions B and D, while 
parametric fitting has been applied to Reduction C. 

The fits we present here are all derived from catalogues 
of sources detected above at least 2.5 a. Thus, although the 
number count estimates are statistical they are fundamen- 
tally different from so-called P{D) analyses in which the 
distribution of pixel flux densities is flt directly to a source 
count model. If we reduced our catalogue to lower and lower 
thresholds we would effectively recover the P{D) results. 
However, this is left to a future paper. 

Tabled] highlights the key steps and differences in each 
reduction's number count estimation. 

5.1 Direct estimate of the differential source 
counts 

The direct estimate works with Reduction D's source list and 
calculates the differential counts directly using the posterior 
flux dens ity distributions fo r individual sources in the list 
following ICoppin et al ](|2005). 

The source list used to calculate the number counts is 
constructed by identifying all 2.5 a peaks in the map, and 
keeping all of the peaks likely to be real (i.e. having < 5 per 
cent deboosted probability of having Si < 0). For the pur- 
poses of measuring the counts the deboosting criterion could 
be relaxed, but with the added complication of statistically 
taking into account the FDR in the counts (see Section [5.21 . 
For this reason, the same criterion that is used to construct 
the SHADES catalogue is applied - the difference being that 
here only Reduction D's data are used. 

An 'effective area' is calculated. It is the area times the 
completeness, and these are estimated together as a function 
of intrinsic source flux density, fl{Si), in order to correct the 
source counts for incompleteness. Fake sources of known flux 
density are injected one at a time into the real maps (with- 
out worrying if the sources fall entirely within the region 
that was measured) and then they are extracted using the 
source extraction method. This procedure is repeated 2000 
times each at flux density levels of 4, 6, 8, 10, 12, 14, 16, 
40, 60, and 100 mJy. A source is considered recovered if it 
is found within 7.5 arcsec of its input position and survived 
the flux density deboosting. Q{Si) is the ratio of the num- 
ber of sources found to the number put in per square degree 
(see Fig. IIH . A smooth best-fitting function of the form 
{S°')/{b + cS") is fitted to the data points and is used in 
correcting the raw source counts. Here S is the fiux density 
and a, b, c axe constants. 

The values of N{Si) in the bins are estimated using 
Monte Carlos, which are also used to estimate the errors. 
In the past, submillimetre survey groups have placed error 
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Table 4. Methods of accounting for bias. Here subscript 'i' refers to the input or unbiased flux density, subscript 'o' refers to observed 
quantities in thresholded maps, and subscript 'd' refers to a true detection. 



Step 


Reduction B 


Reduction C 


Reduction D 


Meaning of poste- 
rior flux density, 
p(5i|5o, (Jo) 


Flux density probability for best fit 
of flux density to the underlying 
'zero-footprint' map distribution. 


Flux density probability for bright- 
est individual source in a measure- 
ment aperture. 


Total flux density in a measure- 
ment aperture. 


Posterior flux density 
expression 


p(5i)p(So,<To|S'i) 


piSi)PdiSo\Si,a^)pi(fJo\Si) 


p{Si)piSa,'ro\Si) 




p(5o,0-o) 




p{So, Co) 


Prior information used 


Actual S/N map for the given field. 


Simulated noisy flux density maps, 
assuming a form of the number 
counts, and completeness for given 
field. 


Simulated noiseless fiux density 
map, assuming a form of the 
number counts, and Gaussian 
photometric errors. 


Source list selection 
criteria 


So/o-o > 3.5 


5o/<To > 2.5 


5o/ao > 2.5, p{Si < 0|5o,(7o) < 
5 per cent 


Number counts 


Place sources in bins at peak pos- 
terior probability. 


Fit prior p{S;) by comparing mod- 
elled p{So,o'o) with data. 


Place sources in bins by integrat- 
ing posterior probability. 


Completeness 


Add individual sources to real 
maps, compare input to output 
catalogue. 


Completely simulate maps, com- 
pare input to output catalogue. 


Add individual sources to real 
map, detect nearest peak and 
compare to input position. 


Spurious Detections 


Completely simulate maps, detect 
nearest peak. 


Completely simulate maps, com- 
pare input to output catalogue. 


See source list selection criteria. 


Counts Uncertainties 


Analytic propagation of errors. 


Monte Carlo, using realizations of 
completely simulated data. 


Monte Carlo, using bootstraps of 
real data. 



bars on the counts by simply accounting for the simple 1 a 
Poisson counting errors (the square root of the raw counts 
in each bin, scaled to unit solid angle). The estimated num- 
ber of spurious and confused sources are sometimes added 
in quadratu re to the lower err or bar in the corrected number 
counts fe.g. IScott et aljf2002fl . Since the deboosting proce- 
dure provides a distribution for each Si, rather than just a 
single value, a modified boot-strapping simulation is used to 
estimate the differential source counts and uncertainties in 
bins of width 2mJy. This simultaneously accounts for the 
Poisson error, as now described. First, a number of sources is 
chosen from a Gaussian distribution centred on the number 
of sources included in the real source list, Ntme, with a stan- 
dard deviation a = \/ iVtruo to account for counting errors. 
This number of sources is then randomly selected from the 
actual source list and their probability distributions sam- 
pl ed with replacement ( i.e. boot-strapping; see Section 6.6 
of IWaU fc Jenkins^l2004^ ; once a flux density is determined 
at random (in proportion to the source's posterior probabil- 
ity distribution), one source per effective area is added into 
the appropriate flux density bin. This procedure is repeated 
10,000 times, in order to make well-sampled histograms of 
the count distributions for each bin. These histograms are 
used to estimate the mean counts and the frequentist 68 per 
cent confldence intervals in each flux density bin and are 
given in Table |3 Simultaneously, the linear Pearson covari- 
ance matrix of the bootstraps across the flux density bins 
can be calculated to assess the correlation between bins and 
this can then be used in model fltting procedures; the co- 
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Figure 11. A fit to f2(Si), the effective area (survey area times 
completeness) of the 850 fim source recovery at each level of input 
flux density for 3 cr map detections which have P{Si < mJy) < 5 
per cent, as determined from Monte Carlo simulations of individ- 
ual sources added to the SHADES maps of Reduction D (see 
Section Ol. The curves are fits to the form (S"')/(6 + cS°-). 
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Figure 12. 850 fim cumulative source counts recovered from a 
fake sky populated with a known source count model (explicitly 
for Reduction D). The diamonds and error bars are the recov- 
ered source counts for these simulated data, observed using the 
same scheme as for the real SXDF data. The ove rplotted his- 
togram (solid line) is the actual realisation of the iBorvs et alj 
l|20n,3r source count model that was used (dot-dashed line). This 
same source count model was used in creating the prior in the 
Bayesian flux density deboosting method described in Section l4.4l 

variance matrix is given in Table |S| for the counts of the 
combined LH and SXDF fields in 2mJy-wide bins. The in- 
tegral counts are obtained by directly summing over the dif- 
ferential counts and are tabulated in Table |S| and shown in 
Fig. EH 

5.1.1 Tests of deboosting on source count recovery 

A test for bias in the recovery of the number counts was 
carried out in the fol lowing way. A fake sky populated with 
the source counts of iBorvs et al.1 J2003I) was created. This 
was observed using the actual SXDF observing scheme and 
a map was made in the same way as for the real data, 
while simultaneously injecting random Gaussian noise with 
an RMS similar to the real map (~ 2mJy). Sources were 
then extracted and deboosted according to the prescription 
described in Section 14.31 The recovered cumulative num- 
ber counts (scaled by the effective area) were found to be 
consistent with the input number count realisation in all 
except the lowest flux density bin, where the uncertainty 
in the completeness estimates dominates (see Fig. I12II . We 
have therefore corrected the differential (integral) counts in 
each of our lowest flux density bins by a factor of 1.6 (1.33). 
The input source counts were also recovered when different 
source count models were used as input to the simulated 
skies, while keeping the form of the prior distribution of 
pixel flux densities fixed. We found that the correction fac- 
tor hardly changed when we used different input models. 



We have also attempted to quantify an overall FDR in 
a different manner to that described in Section | 4.3 | A noise- 
less fake sky was populated with the iBorvs et alJ ll2003l) 
source counts, and the region was sampled using the real 
observing scheme and map reduction steps for each field, 
while simultaneously adding Gaussian random noise into the 
timestream to get the same RMS as the real maps. Sources 
we re extracted in the u sual way and then deboosted using 
the lCoDPin et al] l|2005h prescription to create a final refined 
catalogue of flux density deboosted sources. We found that 
more than 90 per cent of sources detected in the simulated 
maps correspond with input sources above the faintest de- 
boosted flux densities of the actual SHADES catalogue. The 
interpretation of the remainder is complicated, particularly 
as one approaches the confusion regime. We are thus con- 
fident that the overall FDR lies below 10 per cent, but de- 
termining the precise FDR from simulations is complicated 
by source confusion, i.e. interpretation of precisely what the 
source means. 



5.1.2 Another direct estimate of the differential source 
counts 

Direct estimation of the source count density was carried 
out independently working from a catalogue derived from 
Reduction B. The main differences in the approach are listed 
below. 

Instead of calculating an effective area, an explicit coverage 
area. A, corresponding to the portion of a given map with 
observed noise below ctq = 10 mJy is used. Source candidates 
are rejected outside of this region. 

Deboosted source posterior flux probability density func- 
tions are obtained using the normalised histogram of S/N 
in the coverage area of a given map, (So/cto) to estimate 
P(Si\So, ap) in EquationUinstea d of the Gaussian distribu- 
tion used in lCoppin et al^ll2005^ . This is a small difference, 
since the noise is very well described by a Gaussian distri- 
bution. However, the posterior probability density functions 
are truncated outside the region (So/2 < Si < 2So). Exam- 
ples of these deboosted flux density distributions are plotted 
in Fig.Hl 

Source completeness is estimated by Monte Carlo tech- 
niques. This is defined to be C(Si), the detection probability 
for a source of actual flux density Si which is located within 
the coverage area. Sources are inserted into a given map and 
are counted as detected if they are found by the source detec- 
tion algorithm, if their recovered location is within 7arcsec 
of the insertion position, and if the recovered flux density 
is within a factor of two of the input flux density. Simu- 
lated sources recovered within 7arcsec of a genuine cata- 
logue source in this process are discarded. C(Si) is the ratio 
of recovered sources to simulated sources, calculated in half 
mJy bins. 

Source reliability, _R(So/fo), is calculated as a function of 
recovered S/N. For each of the six chop-maps, an artificial 
sk y is generated co nsistent with the source count models 
of Iscott et al.l |2002) and random noise is added, consistent 
with the actual noise maps. Source extraction is performed 
just as it is on the actual survey data. Recovered sources 
are identified with input sources if they lie within Tarcsec 
in position and their recovered flux densities are within a 
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factor of two of each other. R is calculated as the ratio of 
the number of identified sources at a given S /N to the total 
number of recovered sources at that S/N. 

A selection catalogue is formed from Reduction B con- 
taining all sources with S/N > 3.5 and with CTo < lOmJy. 
The mode of the posterior flux density probability distri- 
bution is found for each source, S — mode[Pd(5'i|S'o, (To)], 
and the contribution of each source in this catalogue to the 
number of sources per square degree is calculated as 



AN{S) = 



R{S/ao) 
A X C{S)' 



(3) 



The total source number density in a given flux density bin is 
the sum of AN{S) over all members of the catalogue whose 
mode, S, lies in the flux bin. The uncertainty is calculated 
as the quadrature sum of uncertainties estimated for R, A, 
and C, added to the Poisson uncertainties calculated for the 
number of members of the catalogue with S in the flux bin. 



5.2 Parametric model fits to estimate tlie 
differential source counts 

Number counts are calculated from Reduction C by flt- 
ting models to observed source ca talogues. Thi s tech nique 
is similar to the me thods used by iBorvs et al.l i200 3h and 
iLaurent et ai] ^2005^ for the analyses of SCUBA and BOLO- 
CAM data, respectively. Source catalogues were first gen- 
erated for each field by identifying all 2.5 a peaks in the 
maps. The area of the maps analysed are defined as the re- 
gions having a photometric error CTo < 5mJybeam~'^. The 
model is developed by first expressing the probability dis- 
tribution of the source catalogue ^(5*0, o"o), where So and ctq 
are the observed flux densities and photometric uncertain- 
ties respectively, as the sum of the probabilities of the source 
being a real detection of an object with intrinsic flux density 
Si, Pd{Si, So, fo), and being a false detection, Pf(So, (Jo), 



p{So, (To) 



Pd{Si, So,(^o)dSi +Pf(So,o-o) 



(4) 



The subscript d (f ) is shorthand for the conditional probabil- 
ity that the source is a true (false) detection. Also, note that 
p(So,fo) integrates to 1. The joint probability distribution 
for all of the true detections can be further factored: 



PdiSi, So, (To) = Pd(5'o|5'i, CTo)Pd(cro|Si)p(5i) 



(5) 



The scattering function Pd(S'o|Si, cto) is the probability dis- 
tribution of observed flux densities given an intrinsic flux 
density and measurement error, and hence contains informa- 
tion about the flux bias due to source blending. The function 
J3d(cro|5'i) is the differential completeness, since integrating 
over (Jo is C{Si), the probability of detecting a source with 
intrinsic flux density 5*1 (see Figs. I13| and ll4ll . The flnal fac- 
tor p{Si) is the underlying probability distribution of sources 
with intrinsic flux densities Si. 

The total number of sources detected in the catalogue 
is the number density of sources per square degree A*' multi- 
plied by the survey area A. Multiplying each side of Equa- 
tion 01 by AN, and applying the factorisation in Equation |S] 
gives the observed number density of sources in the cata- 
logue A/' as a function of So and CTq : 



AfiSo,cTo) = Ajpd{So\Si,ao)Pd{cTo\Si)NiSi)dSi (6) 
+ANi{So,ao). 

Here AA(5'o,cro) = ANp{So,ao), N{Si) = Np{Si) is the un- 
derlying differential source counts per square degree, and 
N{ = Apf(5'o,(To) is the spurious detection rate per square 
degree. The left hand side of this equation and the area 
A are measured directly from the survey. Both the scatter- 
ing and differential completeness functions are calculated us- 
ing Monte Carlo simulations, leaving the differential source 
counts N{Si) as the only free parameter to be solved for. 
It is trivial to extend this model to a combined source cat- 
alogue A/''^°™(S'o, o"o) from M independent surveys i taken 
with the same instrument, provided the different p'^{ao\Si), 
Af (5*0,0-0), and A^ are known: 

AA'^°"^(5o,ao) = / Pd(5'ol5'i,ao)A(5'i)^ [AVd('Toi5'i)] dSi 

+Y.A'm{So,ao). 

i 

To calculate pd(5'oj5'i, CTo) and Pd(foj5i) mock bolome- 
ter data were generated using realisations of Gaussian noise 
with the same variance as the real data. To these data were 
added the effect of a population of spatially uniformly dis- 
tributed point sources with a flux den sity distribut i on fol - 
lowing the number counts measured bv lBorvs et alJ i2003l) . 
Source catalogues of all 2.5 a peaks were created in the same 
way as the catalogue for the real data. An attempt was then 
made to identify each observed point source with objects in 
the input catalogue within a 6 arcsec radius. If there were 
intrinsic sources associated with the peak, the brightest was 
considered the match, 5*1. In this way the observed flux den- 
sity in an aperture was related to the flux density for the sin- 
gle brightest source that fell within the measurement aper- 
ture; other fainter sources simply contribute to the upward 
flux bias of this one source through blending. To avoid sen- 
sitivity to extremely faint counts that were not sampled by 
the survey (whose effect is highly model dependent) only 
sources with 5*1 > 2mJy were allowed to be matched (the 
survey was found to be approximately 10 per cent complete 
at this level). The survey is therefore defined to have a com- 
pleteness of for 5*1 < 2 mjy with this model. 

The rate of detection of sources 5*1 , and the distribution 
of Si, So and (Jo using 500 simulated maps was used to esti- 
mate pd( Sol 5i, (Jo) (see Fig. and pd(oolSi) (see Fig. in 
bins of width 1 mJy for Si and 0.125 mJy for So and (Jo. The 
bin sizes for So and Oo were chosen so that the probability of 
having more than one source in a bin is small, so that we can 
use simple Poisson statistics. The coarser bin size for Si was 
adopted so that the Monte Carlo simulations would converge 
more quickly, and since Si does not require high resolution 
because of the width of the posterior flux density distribu- 
tion (see Fig. 1^. At flux densities > lOmJy pd(So|Si, Uo) 
did not fully converge after the 500 simulations since the 
number density of sources is so low (^ 1 source per bin). 
At fainter flux densities it was compared with a Gaussian 
model truncated appropriately for the 2.5 cr source selection 
criteria (see Fig. 1141 and was found to be indistinguishable. 
Rather than run a much larger number of Monte Carlo sim- 
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ulations, the scattering function was instead replaced by the 
smooth theoretical model. Peaks in the map with no intrin- 
sic counterparts were considered spurious detections (either 
pure noise or in extremely rare cases blended sources with 
flux densities < 2 mjy) and were used to estimate the false 
detection rate Nf{So, cto) per square degree using the same 
bins. 

To solve for N{Si) in Equation |7| a discrete non- 
parametric binned model was first adopted, N{j). Each flux 
density bin j was chosen to be 2mJy wide (comparable to 
the photometric uncertainty). Rather than assuming a con- 
stant density of sources across the bin, the counts were mod- 
elled by the product of a free scale parameter, aj, with an 
exponential template function similar to the counts spec- 
trum measured in previous surveys: 



(8) 
(9) 



where Sf^^ and •S'^igh are the flux density limits for the jih 
bin. A downhill simplex optimiser was used to solve for the 
ttj by maximizing the joint Poisson likelihood, jC, of observ- 
ing the true number of detected objects in each bin Af{j, k) 
given the expected distribution Af{j, k) produced by the 
model M{j) in Equation |S| where j and k denote bins of 
So and (Jo respectively: 



C 



n 

j.k 



AA(j, k) 



M{j,k)\ 



(10) 



In order to prevent non-physical answers the aj were con- 
strained to be positive. In addition, it was discovered that 
the solutions were highly unstable and adjacent bins would 
frequently oscillate between and very large values. To rem- 
edy this type of problem the fits were further constrained 
such that the N{j) were monotonically decreasing with Si. 

To calculate the uncertainty in the model, maximum 
likelihood solutions were re-calculated for 500 realisations 
of mock data. These data were generated by drawing the 
same number of sources as in the real list from the maxi- 
mum likelihood distribution N{So,(^o) for the real data. In 
Figs, [ini and the error bars represent the frequentist 68 
per cent confidence intervals for the distribution in each bin. 
The integral source count spectrum (and uncertainty) was 
obtained by directly integrating each model N(j). Finally, 
the 500 fits of N[j) allowed us to directly calculate the sam- 
ple covariance matrix (^N{j), N{k)'^. 

Despite the non-negative and monotonically decreasing 
constraints placed on the binned source counts model, the 
binned differential number counts have a significantly larger 
scatter than was observed in the other groups' estimates. 
This behaviour in the model fitting process is probably due 
to the bin size being inappropriately small given the un- 
certainty in the posterior flux density distributions for in- 
dividual sources, and also the fact that less prior informa- 
tion was used as a constraint (see discussion in Section [5. 311 . 
This problem is analagous to the increased noise one ob- 
tains in an astronomical image when trying to deconvolve 
the point spread function with small pixels compared with 
the FWHM. 
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Figure 14. The scattering function Pti(So|<S'i, cto) at a fixed noise 
level (To = 1.94 mJy for Reduction C. Vertical bands in this plot 
are the distribution of observed flux densities that could be mea- 
sured for a fixed intrinsic flux density and photometric error. In 
the absence of bias and with Gaussian photometric uncertainties 
this distribution is simply a Gaussian with standard deviation 
(To and a mean So equal to Si (the dashed line). However, at 
flux densities Si < 10 mJy the scattering function changes shape. 
The decision to impose a 2.5 (t cut in the source list makes it 
impossible to detect a source at lower S/N and is shown by the 
dotted line. In addition, confusion will tend to cause the faintest 
sources to appear brighter than they really are. Since the 850 ^m 
extra-galactic confusion limit is at ~ Imjy, however, this effect 
is negligible. 



Finally, the model fitting procedure was constrained by 
replacing the binned representation of N{Si) in Equation[7| 
with a smooth parametric model following Equation |21 (for 
consistent notation replace dN/dS with N{Si)). As with the 
binned model, maximum likelihood solutions were found for 
the model parameters N',S',a and 13. The parameter co- 
variance matrix was also obtained using 500 sets of data 
generated from Monte Carlo simulations. The 68 per cent 
confidence envelope for these models is clearly smaller than 
the uncertainties of the individual count bins (see Figs. 1161 
andlTTl. 

All of the analysis undertaken for each field separately 
was repeated using Equation |H] to calculate joint fits of the 
differential source counts to both fields simultaneously. 

Reduction C tested for bias in the recovered number 
counts by simulating 10 data sets with a range of reason- 
able source count models. Source catalogues were produced 
from these data in the same manner as for the real data. 
Using the methods of Section 15.21 (with the same fixed esti- 
mates of pd(<S'o|Si, (Jo) and Pd((Jo|Si)) the recovered binned 
counts were in each case consistent with the input counts 
with insignificant systematic bias. The correction factor in 
the lowest bin is consistent with what was found by Reduc- 
tion D's direct estimate of the source counts (cf. 15.1.11 . 
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Figure 13. The completeness for LH as calculated by Reduction C. The left panel shows the differential completeness pd (ctq | Si ) , the 
probability of a source being detected at a noise level Co given an intrinsic flux density Si. The variations in the vertical direction in 
this plot reflect the relative areas of the map that reached a given noise level. The LH map has a small deep portion with a mean noise 
~ 1.2 mJy, and a large shallower region with a mean noise ~ 2mjy, corresponding to the lower and upper horizontal bands, respectively. 
Marginalization over Co (right panel) gives the more typical definition of completeness, the probability that a source is detected as a 
function of instrinsic flux density, pd('S'i)- 
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Table 5. 850 lira SIfADES differential (in 2 mJy wide bins; quoted per mJy) and integral (in 2 mJy bins) source counts of Reduction D. 
The error bars represent the frequentist 68 per cent confidence intervals of the boot-strapped count distribution in each bin (see text). 
Each differential count flux density bin is indicated by the midpoint of each bin weighted by S~"^, whereas the integral count flux density 
bins are indicated by the lower flux density bound of each bin. The counts and errors in each of the lowest flux density bins of the 
differential (integral) estimates have been corrected by a factor of 1.6 (1.33), for the known undercounting in this bin seen in simulations 
(see Scction l5.1.1i . Note that a Gaussian approximation to the error bars becomes invalid for high flux densities. 



5.3 Summary of differences between the three 
methods 

Having described each group's techniques for calculating 
posterior flux density distributions and number counts, it 
is now appropriate to discuss several key differences: the 
amount of information used; the interpretation of posterior 
flux densities; and the problems inherent in calculating the 
FDR and completeness. 

The philosophies adopted by each group to calculate 



posterior flux density distributions vary in a number of sub- 
tle ways. First, what flux density was calculated? In the 
case of Reduction C the method is attempting to estimate 
the flux density of the brightest source in the beam, while 
Reductions B and D calculate the posterior distribution for 
the beam-smoothed flux density map. The former may be 
more desirable for determining source counts, since it di- 
rectly ties measurements to individual objects. In practice, 
correcting for the confusion is a difficult procedure, which 
depends not only on knowledge of the source counts, but 
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Table 6. Covariance matrix for the 850 /im SHADES combined differential source counts. This can be used along with Table 13 to fit 
models to the counts using the figure-of-merit = i'^ ~ ^)'^ {d — m), where d is the data, m is the model, and is the inverse 

of the covariance matrix. 




Observed Flux Density [mJy] 

Figure 15. Reduction C's observed source distribution (red 
stars) compared with the expected distribution N{So, Co) (shaded 
region) in LH for the maximum likelihood differential counts dis- 
tribution N(Si). The diagonal black line shows the detection 
threshold of 2.5 0^0 that was used to construct the source list for 
determining the counts. The dashed contour indicates where the 
real and spurious parts of Equation|7|are equal. Therefore sources 
detected along this contour have equal chances of being true de- 
tections of objects with intrinsic flux densities greater than 2 mJy, 
or spurious detections (i.e. fainter than 2mjy). Similar to the de- 
tection threshold line in Fig.|7|this contour is not parallel to the 
line of constant S/N. Sources are more likely to be true detections 
in the bottom right region of this plot. 

also their spatial clustering jTakeuchi fc Ishii"2004\ How- 
ever, since the bulk of the sources detected in SHADES have 
flux densities well above the extra-galactic confusion limit, 
simply taking the total posterior flux density in a beam as 
in the case of Reductions B and D is acceptable for the 
purpose of calculating source counts. Furthermore, this de- 
boosting technique avoids the possibility of introducing fur- 
ther model-dependent errors into the posterior flux densities 



in the SHADES catalogue; for this reason the selection pro- 
cedure for the SHADES catalogue uses the ICoppin et alJ 
(2005) algorithm. 

How much information is derived completely from the 
data, and how much information is assumed? Each of the re- 
ductions uses prior knowledge of the 850 /im extra-galactic 
source counts measured in previous blank-field and lensing 
cluster surveys to create simulated source catalogues, and 
to simulate maps exhibiting the chop pattern. Reductions 
B and C use such maps (including noise) to Monte Carlo 
the source detection procedure, and to cross-correlate de- 
tected sources with the input catalogue in order to deter- 
mine completeness. Sources in the output catalogue with 
no corresponding sources in the input catalogue are used 
to measure the FDR based on the selection criteria for the 
source list. Reduction D uses the simulated maps (with no 
instrumental noise added) strictly as a prior for the poste- 
rior flux density distributions. Spurious sources are handled 
by constructing the catalogue in such a way that it is nearly 
free of spurious sources. Completeness is determined by in- 
troducing fake sources into the real map and attempting to 
recover them over a range of input flux densities. An ad- 
ditional completeness correction is calculated in a manner 
similar to Reductions B and C by ensuring that the recov- 
ered source density matches the input source density in a 
separate simulation that includes both sources and noise. 

Reductions B and D both choose to calculate the differ- 
ential source counts in bins directly from the posterior flux 
density distributions of sources in their catalogues. There is 
clearly a possiblity for the prior used in the posterior flux 
density calculation to bias the estimate of the source counts. 
Reduction B handles this problem by using the pixel dis- 
tribution in the observed map itself as a prior for individ- 
ual sources. Reduction D tests the procedure on maps with 
sources drawn from several different source count models, 
while keeping the prior for de-boosting flxed to test for bias 
(see Section lS.l.H . Rather than assuming a prior, Reduction 
C leaves the prior as a free parameter and attempts to flt it 
to the observed catalogue. 

The simulations undertaken by all reductions indicate 
that the uncertainty in the completeness introduced by dif- 
ferent source count models is dominated by other uncertain- 
ties in the counting procedure (map noise and Poisson count- 
ing uncertainties). In the case that sources are completely 
isolated, the probability that they are detected (i.e. com- 
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pleteness) is only a function of the noise. On the other hand, 
source blending near the detection threshold may affect com- 
pleteness as a function of the underlying source density. 
Sources that blend together and are detected as a single 
bright source increase the completeness if individually they 
would not have been detected, but decrease the completeness 
if individually they could have been detected. In this paper, 
the faintest sources that we claim to count are > 2 mjy. The 
full range of integral source counts at this flux density mea- 
sured here, and in previous work, is conservatively between 
1000 and 10000 sources per deg^. Therefore, for spatially 
uniformly distributed objects there are on average ~ 0.02- 
0.2 sources that land within a beam (the solid angle of the 
SCUBA beam is ~ 2 x 10"^ deg^; the FWHM^). Given the 
Monte Carlo simulations undertaken in Sections lS.ll and 15.21 
using a plausible range of different input source counts mod- 
els, it is therefore not surprising that the variations in the 
completeness corrections down to 2mJy are dominated by 
uncertainties other than the variations in the counts model. 



Which method is best? Since the 850/im source counts 
are known well enough to construct a useful prior, it makes 
sense to use this information in the analysis of SHADES 
data. It is not surprising that Reductions B and D quote 
smaller uncertainties than Reduction C in the binned dif- 
ferential counts, since more information is being used. On 
the other hand, for surveys at a wavelength for which little 
information is known, the technique for binned counts de- 
scribed in Section 15.21 is more conservative. The technique 
adopted by Reduction C is more useful once a smooth para- 
metric model is assumed to describe the source counts. Con- 
straining the counts in this way produces a range of mod- 
els, with a spread generally consistent with the smaller er- 
ror bars quoted by the other groups for the flux density 
bins. Furthermore it offers the cleanest way for constrain- 
ing the model; the model parameters are varied directly to 
calculate the likelihood of observing the source catalogue. 
With Reduction D, on the other hand, parametric models 
are fit to the binned counts calculated previously. However, 
the binned differential counts from Reduction D can more 
easily be combined with counts from other surveys to con- 
strain models over wider ranges in flux density, because the 
data products that come out of the procedure are binned 
counts and a bin-to-bin covariance matrix. 



For this last reason we proceed with the number counts 
and bin-to- bin covariance matrix of Reduction D. It is rea- 
sonable to select one group's reduction, since we cannot eas- 
ily combine the three sets of counts like we did for flux den- 
sities or positions, and moreover all of the counts across re- 
ductions appear to be consistent with each other within the 
error bar estimation, for all but the lowest bin (where Reduc- 
tion B appears to be low). In model tests we quote the best 
fit parameters using the counts of Reduction D, since that 
reduction provides the smallest error bars (cf. Section f4. It . 
while using the 68 per cent confidence intervals of Reduction 
C as a consistency check. The differential counts for LH and 
SXDF are shown in Figs. I16i and ll7l respectively. 
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Figure 16. Differential souree count densities in the LH region. 
The number of sources deg~^ mjy~^ at a given flux density is 
plotted against flux density (with different reductions offset hori- 
zontally for clarity) . All error bars are estimates of 1 cr uncertain- 
ties of the full error distribution including completeness estimates 
(see text). The horizontal line marked 'Survey Limit' corresponds 
to the 68 per cent confidence limit which can be drawn from find- 
ing zero sources in 1/8 dog^ (in a 2mjy-wide bin). The 68 per 
cent confidence interval on the set of best-fitting double power- 
laws of the form of Equation |^ (fit to the differential counts of 
Reduction C) are overplotted as dashed lines; the error bars show 
that there is little constraint above this level. 



6 RESULTS: MODELS AND CUMULATIVE 
SOURCE COUNTS 

6.1 Fits to differential counts 

In fitting models of the numbers of sources at different fiux 
densities to our data there are different approaches which 
could be taken. The first choice is whether to flt to the list 
of sources directly or to count the sources in a series of dis- 
crete flux density bins and flt to those numbers. In principle 
these methods could each be applied to the results of each 
reduction, but there is little to be gained from generating 
eight different flts in this way. We have tried one approach 
each on the data of two separate reductions. 

Starting with the catalogue from Reduction C, we 
use Monte Carlo techniques to compute the maximum- 
likelihood distribution of models that fit the source cata- 
logue, given the completeness and FDR. The dashed curves 
in Figs. ITHl and 1171 indicate the 68 per cent frequentist in- 
terval on the best-fitting models. This procedure requires 
detailed knowledge of the FDR and completeness as a func- 
tion of flux density and therefore it is difficult to compare 
directly with published counts. Note that the 68 per cent 
range on the models flt to Reduction C is about the size of 
the individual error bars in Reduction D, and that the re- 
ductions are all consistent with each other within the error 
bars in all but the lowest flux density bin. 
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Figure 17. Differential source count densities in the SXDF re- 
gion, as for Fig. 1161 Notice that in all except the lowest two flux 
density bins for Reduction C, the density of sources inferred in 
the SXDF is approximately 1/4 0- lower than for the LH data; see 
Section l6.6l for a discussion of possible field-to-field variations. 

Starting with the counts and covariance matrix, based 
on Reduction D (found in Tables |^ and (SJ , we fit a simple 
power-law of the form 



using a minimum parametric approach. Here, A'^' is the 
normalisation at S' = 5mJy (the approximate 'pivot' point 
for our data). 

The 'Survey Limit' (as indicated by a dot-dashed line in 
the figures) corresponds to the 68 per cent confidence limit 
which can be drawn from finding zero sources in 1/8 deg^ (in 
a 2 mjy-wide bin) in each of the two fields. This provides a 
rule-of-thumb counting limit to which the SHADES fields 
are sensitive, given the observed areas and depths achieved. 
Note that the Gaussian error approximation breaks down for 
the highest flux density bins but these carry little weight in 
the fits in any case. The approach taken by Reduction C also 
makes full use of the range of fiux densities in which there are 
no detected objects, and thus provides a consistency check. 

The best fitting line has parameters: a — 2.9 ± 0.2 and 
N' = 189 ± 26 in LH; and q = 3.0 ± 0.3 and N' = 136 ± 24 
in SXDF. Note that a and A''' are essentially uncorrelated 
with each other, because we choose the normalisation to be 
around S' = 5mJy, near the centroid of our data. This re- 
sult, that a is virtually the same in both fields while A'^' is 
lower in SXDF, is consistent with the data in Figs. ITCl and 
1171 where every bin in SXDF appears lower than the corre- 
sponding bin in LH by about 1/4 cr uniformly across source 
brightness. The slope result agrees with previous estimates 
obtained by other grou ps (we quote 1 a limits from the lit- 
erature): a = 2.8 ± 0.7 lBlain etHl JTgQQ.') : a = 2.9 ± 0.25 



Figure 18. Differential source counts in the combined fields (LH 
+ SXDF). The number of sources mjy~^ deg~^ at a given flux 
density is plotted against flux density. All error bars are estimates 
of 1 cr uncertainties of the full error distribution, including com- 
pleteness estimates (see text). The horizontal line marked 'Sur- 
vey Limit' corresponds to the 68 per cent confidence limit which 
can be drawn from finding zero sources in 1/4 deg'^ (in a 2mjy- 
wide bin). The best-fitting single power-law fEguation llll dashed 
line), Schechter function fEauation 1131 dotted curve), and bro- 
ken power-law ('Eauation ll2l solid lines) to the differential counts 
of Reduction D are plotted. It is clear that the Schechter func- 
tion or the broken power-law are a better fit to the data than 
the single power-law, indicating a break somewhere in the range 
~ 5-13 mJy. 

iBorvs et alj tOO^ : a = 3 2+5 ^ iBarger. Cowie fc Sanders! 
^1999^ : and a ^ 3.25 ± 0.7 lEales et alJ tOOd) . 

The total number of sources per square degree in LH 
versus SXDF differ by less than 2 a (for purely Poisson scat- 
ter) and there is no convincing evidence to suggest that the 
ratio of the counts between the fields is different from 1. 
In the rest of our analysis we combine the data for the two 
fields to improve the statistical power slightly. The combined 
counts (i.e. for both fields) are plotted in Fig. 1181 

It is clear that a single power-law model for the com- 
bined data, shown as the dashed line in Fig. 1181 is a poor 
fit to the data. The fit remains poor even if the counts in 
the lowest bin are arbitrarily doubled, and we do not believe 
we have mis-estimated completeness by this much. We have 
thus explored fits to broken power-laws, and to a Schechter 
function ( Schec hter 1976 ) . 

A full fit to Equation |5| involves solving for four pa- 
rameters: A'^', S' , a, and fH, which are not simultaneously 
well-constrained by the few data points. The parameter er- 
ror bars are therefore large and correlated (see Table [7J. If 
instead, we hold the break flux density fixed at S' = 9, a 
visually plausible value (and close to the best-fitting value of 
S'), the resulting uncertainties in A'^', a, and (5 are reduced, 
as also shown in Table |7| and the errors become almost un- 
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Equation 


S' 


N' 


a 


/3 




9.4 ± 4.0 


647 ± 739 


2.0 ± 0.6 


6.0 ± 2.2 


m 


9 


735 ± 123 


2.0 ± 0.2 


5.8 ±0.9 




9 


49 ± 9.3 


2.4 ± 0.2 


5.1 ± 0.9 




3.3 ± 1.2 


1599 ± 1183 


-2.0 ± 0.7 





Table 7. Best-fitting parameters of Equations |3] (smooth double 
power-law). 1121 (a broken power-law) and 1131 fSchechter function) 
fit to the combined counts. Fixing the break location in the source 
counts, 5", has the effect of reducing the error bars on the other 
parameters (cf. top 2 rows). 



correlated. The of the broken power-law fit is reduced 
from the single power-law by about 10 for one new parame- 
ter (A'^', Q — > A^', a, /3). However, we note that the resulting 
reduced x values are unrealistically small (i.e. less than 
1.0) in all but the single power-law fits. This probably re- 
flects the fact that the errors are non-Gaussian or that there 
are some correlations not taken into account properly in the 
covariance matrix. But the fact remains that fits to a single 
power-law do not describe the data well. 

We also fit a broken power-law of the form 



dN 
dS 



N' 



Jt) forS > 5" and 



(12) 



forS < S', 



using the same minimum parametric approach. The best- 
fitting form of Equation 1121 holding S' fixed at 9 mjy, is 

plotted in Fig. llSl as 2 sol id lines. 

A Schechter function iSchechterll97d) was also fit to the 
counts, since it is at least a physically- motivated functional 
form: 



dN N' f S . , 



(13) 



The parameters TV', S' and a of the best-fitting form of 
Equation 1131 are tabulated in Table |7| (see Fig. 1181 dotted 
curve) . The Schechter function above fits the data as well as 
the other 2 power-laws we fit here. This may not be surpris- 
ing given that the relationship between 850 /im flux density 
and rest-frame FIR luminosity is nearly constant and inde - 
pendent of the redshift for 1 < 2 < 8 (see lBlain et al.ll2002l) . 
But there is no evidence to favour an exponential fall ofi^ 
over a steeper power-law at large flux densities. 

The data clearly favour a change in the slope for the dif- 
ferential counts, although the precise position of this break 
is not well constrained. This change in slope of the differ- 
ential counts should be helpful in breaking degeneracies in 
fltting models of luminosity function evolution. 

If we use additional information from low flux density 
counts (particularly from cluster lens fields) the evidence for 
a break becomes stronger. Also, a shallower slope at low flux 
densities is required in order not to overproduce the suhmil- 
limetre background. However, neither of these arguments 
requires that the break be at a large enough flux density 
to be seen within the SHADES data. It is our direct fits to 
the differential counts which, for the first time, constrain the 
power-law break to be at a flux density of several mJy. 



6.2 Background estimate 

One can use any of the best-fitting models to estimate 
the total flux density which has been resolved into dis- 
crete sources. Adopting Equation 1121 and using the best- 
fitting parameters in Table |7| with 5*' fixed at 9 mJy, we 
can calculate the total 850 flux density by integrating 
S N{S) dS . The estimated background of sources brighter 
than 2mJy is 9.7l3;^ x lO^mJydeg-^ We can also esti- 
mate the background directly from the SHADES counts by 
performing the sum SN{S) over the bins, which gives 
l.Otoil X 10''mJydeg-2. 

By comparison, the total FIR 850 back- 
ground inferred fro m COBE - FIRA S is 3.1- 
4.4 X lO'^mJydeg"^ dPueet et al.l Il996l : iFixsen et all 
119981; iLaeache. Pueet fc Dole '20051 . We flnd that a survey 
complete down to 2 mJy at 850 ^m would resolve between 20 
and 30 per cent of the FIR background into point sources. 
The uncertainty in this fraction is dominated by uncertainty 
in the unresolved background, rather than uncertainty in 
the SHADES sources. Our result is consistent with values 
quoted by other groups (e.g.lBarger. Cowie fc Sandersll999l . 
lEales et alJ bood iBorvs et alJl2003^ . A significant fraction 
of the suhmillimetre emission lies below the detection 
limit of blank-field SCUBA surveys at 850 fj,m. Therefore, 
knowing the number counts accurately down to much 
fainter flux density limits (~ 0.1 mJy) is essential in order 
to constrain models that predict the evolution of luminous 
IR galaxies. 



6.3 Cumulative source counts 

Previous SCUBA surveys with fewer sources have reported 
source densities in the form of cumulative counts rather than 
differential counts, since the latter are somewhat sensitive 
to the bin choice for small source catalogues. The SHADES 
integral counts are obtained by directly summing over the 
differential counts of Reduction D. The best-fitting double 
power-law model and Schechter function (cf. Equations 1121 
and 1131 and Table |7J of the differential counts have been 
integrated to produce models of the cumulative counts. The 
cumulative counts and these models are plotted in Fig. 1191 
and given in Table |5] Note that the models are not fit in 
the cumulative counts domain. These should be the most 
accurate 850 /im number counts in the flux density range 
2-15 mJy achieved so far by any single survey. 



6.4 Comparison of cumulative counts to previous 
estimates 

Previous measurements of the cumulative number counts 
have differed by factors of ~ 2.5-3 among groups ob- 
serving various small (10s to ~ 200arcmin 
of sky, especially at flux d e nsiti es between 

_ 199' 

200d . IScott et alJ KQO: 



1 



e .g. iBarg er. Cowie fc Sander^ 
IJ l200d . IScott et alJ 



■^^ patches 
2- 6 mJy 



iBlain et al 



Eales et alJ 



199£ 



200s 



Webb et al J 120031: see discussion in IScott et al.l 120021) . The 



culprits are most likely sampling variance, clustering and/or 
low number statistics, due to the small-area surveys, as well 
as different ways of treating flux density boosting in the low 
S/N regime. 

We might expect most agreement at the bright end 
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Figure 19. Cumulative combined SHADES source counts (solid points) compared to previous estimates. The 68 per cent confidence 
interval in Reduction C's model fits is indicated by the dashed lines. The 95 per cent Poisson upper confidence limit to the surface 
density of sources brighter than 22 mJy in SHADES from Section 16.51 is shown by the downward arrow. The source counts deter- 
mined from various other blank -field and cluster lensing surveys (corrected for lensing) are indicated using different symbols, with the 
IScott. DunloD fc Seriean3 i2006 'l compiled counts being represented by the hatched region. The best fitting form of Equation 1121 fdark 
solid curve), holding S' fixed at 9mjy, and Equation 1131 (dark dotted curve) to the differential counts are integrated and plotted here. 
See Table 171 for the best-fitting parameter values for these functions. The data confirm a break in the power-law somewhere in the middle 
of the flux density range. The lensing data are consistently high, which could have a number of explanations (see text). 



of the source counts (among previous wide-area surveys to 
varying depths), since only the brightest sources would have 
prevailed. However, we might also expect most disagreement 
between small area blank-field surveys, since very bright 
sources are rare and a large degree of sampling variance may 
therefore come into play. In Fig. 1191 the SHADES counts are 



plotted in comparison with some previous work in order to 
assess the level of agreement. 

Gravitational lensing by clusters of galaxi es has been 
used as a tool to study the fain t est SMGs (e.g. |Blai net_alJ 
[ 19991: ICowie. Barger fc KneibI 120021: ISmail et al.l l2002l: 
IChapman et alJl2002l: IWebb et alJl2005l: and lKnudsen et alJ 
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l2006t) . The agreement between these independent surveys 
is very good below about 2mJy. We note that there is a 
smooth transition also at 3mJy betwee n the cluster s urvey s 
and this wo rk. However, we agree with lBorvs et al] i2003l) . 
IWebb et alJ JjoOSn and others who noted that the lensing 
number counts above 4mJy (where the cluster and blank- 
field counts overlap) appear higher than the combined blank- 
field survey counts. This is most likely due to the fact that 
the quoted error bars in the literature do not contain a vari- 
ance component and that the flux density boosting bias has 
been treated differently by different groups. 

At the bright end, the iBorvs et all J2003l) blank-field 
counts follow an approximately power-law decline with in- 
creasing flux density, whereas ours appear to steepen more 
dramatically, but t he results agree we ll at low flux densi- 
ties (2-6 mjy). The lBorvs et al.N2003t) counts appear to be 
higher on average than the SHADES counts and those of 
IScott. iJunloD fc Serjeant (2006), but are within the error 
bars. Sampling variance may be responsible for these differ- 
ences, since no survey so far has enough area to overcome 
this, particularly if these objects cluster on arcminute scales. 
SHADES agrees well, pa rticularly at the high er flux den- 
sities, with the counts of IScott. Dunlop fc Serje ant ( 2 00^ . 
which is a compilation of data covering a similar area to 
each of the two SHADES flelds (the 'band' in Fig. [T^ . 
For flux densities > 10 mJy the S HADES counts ap pear 
to be inte rmediate between those of lBorvs et alJ i2003ft and 
IScott. DunloD fc Serjeant, 1.2006i') . 

6.5 Bright source constraint 

We can also estimate a limit to the surface density of the 
brightest SMGs. There are various ways to do this, but 
the simplest is just to take the fact that SHADES con- 
tains no sources brighter than 22 mJy in the entire surveyed 
area to constrain the bright counts using Poisson statis- 
tics. We find A''(>22mJy) < 17deg"^ at 95 per cent con- 
fidence. This estimate is probably more robust than previ- 
ous upper limits at sim ilar flux densities (e.g. Borvs ct al. 
l2003l: IScott et alJ |200^ ari d can be compared with the 
even brighter iBarnard et all ||2004) limit of iV(>100 mJy) < 
2.9 deg-^ 

6.6 Comparison of the two fields — evidence for 
sampling variance? 

We can estimate the expected sampling variance between 
the flelds in the following way. The va riance of counts in 
cells is generally given by jPeebleslll98nl. eqn. 45.6) 

(j^{M) ^ MVi + j j vj{ei2)dQ.^dQ.2, (14) 

where N is the number per unit area on the sky, Q. is the 
solid angle of the cell, ■w{9) is the angular 2-point correla- 
tion function and the integrals are over all separations 6\2 
between 2 positions. If the size of the flelds is much bigger 
than the correlation length, then one of the integrals can be 
done trivially, giving 

o'^{N)^Nn{l+N jw{e)dQ.], (15) 

where 6 can now be considered as the angle from the cen- 
tre, say. The term with the integral is the excess variance 



over Poisson. If we are in the power-law regime for w{6), 
then wiO) = {9/9o)~'' , with 7 ~ 0.8 typically. For simplicity 
we can assume that the area is a circle of radius Smax. In 
that case the integral gives (27r&o */1.2)^i;ax (with any cut- 
off at some Omin being relatively unimportant). Eventually 
this will break down, since w{6) goes negative, but it is a 
reasonable approximation in the power-law regime. Putting 
in n ~ 400 arcmin^ and AfQ ~ 60 for SHADES, one obtains 
an RMS scatter in the counts of 

a{jV) ~ {JVny^^ X (1 + 0.5{eo/l"fy^^. (16) 

Measuring the value of the clustering angle, Oq is one 
of the key goals of SHADES, and a careful analysis will be 
made once the forthcoming redshift distributions (Itziar et 
al. and Clements et al., in preparation) are included. One 
might ask what constraints on do can be obtained from com- 
paring the number of sources per square degree found in the 
two SHADES flelds. In SectionlOwe found 189± 26 sources 
deg~^ mjy~^ at 5* = 5mJy in the LH and 136 ± 24 in the 
SXDF. The uncertainties are essentially Poisson, and are 
therefore underestimated by a factor 

J^{eQ) = {l + 0.5{9o/l")°y^^ 

if clustering is important. In Fig. 1201 we have plotted the 
likelihood to observe a difference, AA*' between the num- 
ber counts in the two fields which is larger than the one 
we have observed as a function of clustering angle, 9o. 
Ivan Kampen et al.l (|200^ estimate that for surveys like 
SHADES, 9q might range from 5arcsec (for a hydrodynamic 
model) to perhaps 20 arcsec (for a high mass merger model) . 
Notice that the likelihoods in Fig. l2Ul remain within the cen- 
tral 68 percent region even over a broader range of 9o than is 
anticipated, corresponding to values of .7-" running from 1 to 
~ 4. Thus, for this particular comparison we do not provide 
even an interesting 1 a constraint on the clustering angle. 



7 CONCLUSIONS 

This is the first instance that a SCUBA data set has been 
subjected to a detailed comparison of differing independent 
reductions. As a result, we believe we have produced the 
most robust and reliable submillimetre-selected source cat- 
alogue to date. We have learned that using monthly average 
FCFs is just as good as using nightly calibration observa- 
tions, but with the advantage of less time lost to calibration 
observations. Also, using larger pixels (3 arcsec as opposed 
to 1 arcsec pixels) to make maps of these data does not seem 
to increase the positional uncertainty in the centroiding of 
detected sources. 

The SHADES survey is the final legacy of SCUBA, be- 
ing the largest extragalactic survey undertaken to investi- 
gate the nature and redshift distribution of a complete and 
consistent sample of > 100 SMGs. Here, we have presented 
850 nm source catalogues and number counts (differential 
and integral) using information from four independent re- 
ductions of the data, complete to about 720arcmin'^, down 
to an RMS noise of ~ 2.2 mJy. We also provide 3 a up- 
per limits to the 450 fj,m flux densities of each SHADES 
source. With approximately double the area of all blank- 
fleld surveys combined, this work has produced the most 
robust suhmillimetre source catalogue for multi-wavelength 
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Figure 20. The likelihood that fluctuations would lead to a larger 
difference in the density of 5 mjy sources than has been observed 
between the LH and SXDF regions is plotted against clustering 
angle, Oq. Even though the run of Oq which is plotted is broader 
than the range anticipated from physical models, the likelihood 
values all lie in the central 68 per cent of the distribution. Mild 
clustering is a slightly better fit to the data than no clustering 
{9o = 0), but these data do not provide an interesting constraint 
on clustering amplitude, even at the 1 a level. 



estimating redshifts and determining luminosities and other 
physical properties for the sample. A detailed analysis and 
modelling of the number counts will be presented in fu- 
ture paper (Rowan- Robinson et al. in preparation) . Together 
these will provide much tighter constraints on galaxy evo- 
lution models. We note that wider surveys, such as those 
planned for SCUBA-2, the LMT, Planck or Herschel, wiU 
be able to place better constraints on the bright end of the 
number counts, since these SHADES maps are insufficiently 
wide to detect rare bright SMGs. On the other hand, deeper 
surveys, such as those also planned for SCUBA-2, are needed 
to push into the JCMT's confusion regime, and eventually 
ALMA will allow us to study these SMGs with the desired 
angular resolution. 
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studies, as well as the most accurate measurement of SMG 
number counts. 

We find that the differential counts are better fit with 
a broken power-law than a single power-law, the location of 
the break in the source counts being at several mJy. Using 
the ratio of the counts between the two SHADES fields, we 
are unable to provide a useful constraint on SMG clustering. 

Overall, the increased accuracy of the source counts 
due to the new measurement by SHADES will allow 
us to study the evolutionary nature of SMGs in much 
more detail than has previously been possible. SMG num- 
ber counts have been able to put severe constraints on 
galaxy evolution models, where new ingredients or as- 
sumptions had to be adopted in order to reproduce them 
(e.g. Kaviani. Hachnclt & Kauffmann 2003, Granato 2 004 . 
feaugh et akfcooa) . 

The maps and catalogues that have resulted from this 
work, in combination with multi- wavelength measurements 
from deep radio, far-, mid-, and near-infrared, optical and 
X-ray imaging, will yield the photometric redshifts neces- 
sary to investigate the SMG contribution to the cosmic star- 
formation history of the Universe and to discover the true 
nature of the SMG population. 

Further work on the SHADES data wiU include a P{D) 
fluctuation analysis (which should be relatively straightfor- 
ward because of the fairly uniform noise), and clustering 
estimates (again, significantly easier for SHADES than for 
previous surveys), as well as follow-up studies, focussed on 
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APPENDIX A: NOISE SPIKE REMOVAL 
TESTS 

In order to assess its effect on the data, an aggressive at- 
tempt was made to remove the noise spike (see Section]^ as 
follows. First, the array average was removed as a function 
of time, using only the least-affected bolometers to estimate 
the sky at each time, so as not to introduce the spike into 
bolometers which do not already exhibit the excess power. 
Next, another function of time (essentially a template of the 
spike effect in the timestreams) was subtracted, for which 
each bolometer sees a different fraction. This fraction was 
found by minimising the spread in each bolometer over time 
(i.e. using a minimum variance estimator) and is expected to 
be close to for the bolometers which exhibit little evidence 
of any excess power and close to 1 (after being appropriately 
normalised) for the bolometers exhibiting the highest degree 
of excess power. This approach was applied to Reduction D, 
but the results should be applicable to all reductions. 

Approximately 20 per cent of the SXDF data files ap- 
pear to be contaminated by the noise spike. While all 37 
bolometers may contain the effect to some extent, about 
30 per cent of them exhibit a significant excess of power in 
the power spectrum at a period of 16 samples. So overall, 
about 8 per cent of the data going into the SXDF map is 
significantly affected by the noise spike. As a result of the 
spike removal, the RMS of Reduction D's SXDF map was 
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decreased by a mere 1 per cent, as compared to the reduction 
in which the spike's presence was ignored. 

We find that the LH map is shghtly more affected than 
the SXDF map: 30 per cent of the Lockman data files are 
affected and once the spike is removed the map RMS is im- 
proved by 3 per cent. For the 10/14 Reduction D-detected 
LH SHADES catalogue sources lying in regions with more 
than 50 per cent of their data affected by the spike, no po- 
sitional offsets are detected from the corresponding spike- 
removed map sources. Positional offsets of up to the 3arcsec 
pixel size occur in the spike-removed Reduction D map for 
only 4 of the LH SHADES catalogue sources. 

In general the contaminated files are distributed fairly 
uniformly around the maps, although the worst examples 
are localised to specific regions. Hence even although the 
overall map RMS changed very little, it may still be that 
the flux densities or positions of particular sources could 
be affected significantly. We checked for signal, noise and 
S/N variations and found that small differences do occur 
for sources detected in the noise-spike removed maps versus 
the regular maps, as expected. The flux densities and S/N 
values differ from the regular map values by about 5 per cent 
on average, and never by more than 10 per cent. Overall 
28 per cent of the SHADES sources have been affected to 
any discernible extent by the presence of the noise spike, 
in terms of flux density, S/N or position. Approximately 25 
per cent of these sources show positional offsets of < 3- 
4 arcsec from their original positions in the map, while the 
rest show no discernable positional differences. Therefore, we 
conclude that the noise spike mainly contributes additional 
random noise to the source flux densities and centroids of the 
data, and thus can be safely ignored in our data reduction 
treatment. However, for a few of the most-severely affected 
sources we have included a discussion of the effects in the 
individual source notes in Appendix lUl 



APPENDIX B: TESTS OF THE DATA 

Several different checks were performed by the reduction 
groups in order to examine the noise properties of the maps 
and test the individual source robustness by examining spa- 
tial and temporal fits to point sources in the maps. 



Bl Tests of Gaussianity in the maps 

Reduction C's photometric error distributions were tested 
for Gaussianity by performing Monte Carlo simulations. A 
realisation of Gaussian noise was created from the measured 
detector signal variances for SXDF. The simulated noise sig- 
nals were reduced and rebinned into pixels using the same 
procedures as for the real data. Sources were placed into the 
map one at a time at random locations and their recovery 
was attempted. The joint probability distribution from the 
simulations was calculated using the same technique that 
was used for the real data. We found that the photometric 
uncertainty is very close to Gaussian, as expected, and the 
parameter, a that comes out of the source extraction routine 
is entirely consistent with what is expected. 



B2 General source robustness 

A quick test of source reality is to create the negative of 
the map and search for sources using the same triple-beam 
template. Aside from pixels associated with the off-beams of 
positive detections, we find 3 and 4 > 3.5 a 'detections' and 
21 and 25 > 3.0 a 'detections' in the inverted reduction D 
SXDF and LH maps, respectively. This is consistent with the 
expected number of false positive detections in noisy data. 
Note that only one of these negative sources (3.9 a, in LH) 
has high enough S/N and/or low enough noise to survive 
the deboosting procedure of Section f4.4l 

We have also performed spatial and temporal tests in 
order to determine how well the raw timestream data fit the 
final PSF-fitted maps. A candidate source may be 'detected' 
in the map, but it may not necessarily be well-fit by the 
PSF, or it may be a poor fit to t he set of differ e nce d ata, 
or both. See lPope etaP (l2005l) and lCoppin et al.l (l2005l) for 
more discussion of such tests in other SCUBA surveys. Note 
that these two tests also check for Gaussianity of the noise 
in the maps, but that they are not a strong test of this 
distribution. 

The 'spatial x^' test provides a gauge across the map of 
the goodness-of-fit of the triple-beam differential PSF to the 
data, and thus indicates if a source is poorly fit by the as- 
sumed PSF. We have picked out > 2.5 a sources with val- 
ues lying outside of the ±2cr = 2.16 (LH), or 2.05 (SXDF) re- 
gions of the complete spatial distribution of the smoothed 
LH and SXDF maps. We find 6 such sources in Reduction 
D's LH map, 3 of which correspond to SHADES catalogue 
sources, and we find 16 such sources in Reduction D's SXDF 
map, 5 of which correspond to SHADES catalogue sources. 

The 'temporal x^' provides a measure of the self- 
consistency of the raw timestream data which contribute 
flux density to each map pixel. We have selected > 2.5 cr 
sources in the maps lying outside the ±2(j= 1.60 (LH and 
SXDF) regions of the temporal x^ distribution. In Reduc- 
tion D's LH map, we find 17 such sources, 5 of which cor- 
respond to SHADES catalogue sources. In Reduction D's 
SXDF maps, we find 16 such sources, 6 of which correspond 
to SHADES catalogue sources. 

We have noted those source candidates with relatively 
poor spatial and/or temporal x'^ values in Appendix IHl All 
of these flagged sources have either temporal or spatial x'^ 
values just outside the 2 a regions of the x^ distributions, 
except for one source that has a spatial x^ value > 3 cr 
(LOCK850.15). 



B3 Analysis of the 450 /^m data 

450 fim data were taken simultaneously with those at 
850 fim, providing complementary short-wavelength pho- 
tometry (or flux density limits) over a poorly sampled wave- 
length regime of the spectral energy distribution of SMGs. 
For objects in the SHADES catalogue, these data can pro- 
vide useful constraints for Far-IR-to-mm wavelength photo- 
metric redshift estimates, a powerful tool for measuring the 
star formation history of the submillimetre galaxy popu- 
lation (e.g. Aretxaga et al. 2003). However, 450 /im obser- 
vations of faint objects with SCUBA are of limited use, 
as: (1) the JCMT beamshape is non-Gaussian at 450 /^m, 
due to the non-optimization of the telescope surface for 
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short-wavelength observations; (2) the atmosphere is more 
opaque at this wavelength and therefore the noise will be 
more sensitive to small variations in the sky transmission 
than at longer wavelengths; and (3) atmospheric emission 
fluctuations are much more severe at 450 /im, and are not 
sufiiciently reduced by removing the 1 Hz array average 
(which is so successful at 850 /xm). As the SHADES survey 
was conducted throughout a range of weather conditions 
(tcso > 0.05; mostly unsuitable for sensitive 450 /^m obser- 
vations), these data are not expected to add much to our 
understanding of the 850 /xm source sample. Nevertheless, 
for completeness we describe the main results derived from 
our analysis of these data. 

B3. 1 450 /im reduction strategies 

Reduction methods applied to the 450 /im data are similar 
to those described in Section |21 with the following minor 
exceptions. Reduction C uses smaller, 1 arcsec pixels when 
rebinning, and adopts a 7.5 arcsec Gaussian rather than the 
full, chopped PSF to extract flux densities. Reduction D 
also ignores the off-beams in the data and thus does not 
fold them in when rebinning. The results of the Reduction 
A analysis are not included in the 450 comparison. 

The spatial variation in 450 /im flux densities across the 
rebinned maps is more correlated between Reductions B and 
C than between either of these and Reduction D. This is not 
surprising, given that Reductions B and C follow a similar 
calibration strategy, i.e. applying FCFs derived from cali- 
bration observations taken on the same night as the data, 
rather than adopting a monthly average FCF as was the 
strategy for Reduction D (cf. Table 0. There are plausible 
reasons why using frequent measurements of the FCF can 
either be beneficial, or detrimental, and our data cannot dis- 
tinguish between these possibilities, as we do not know the 
true amplitude of the underlying signal. Based on the con- 
sistently lower noise in reduction B, one might conclude that 
the longer, monthly average calibration strategy of Reduc- 
tion D is introducing noise, and that at 450 ^m the more 
frequent calibration strategy adopted by Reductions B and 
C is optimal. However, within the large uncertainties of the 
450 pim photometry (described below), the reductions are 
broadly consistent. 

B3.2 A search for 450 blank-field sources 

As the beamsize at 450^m is roughly half that at 850 ^m, we 
expect to uncover many more spurious sources in each map 
above a given S/N threshold. This translates to a greater 
number of blank-fleld source candidates (i.e. those found 
at random positions in the map), so care must be taken to 
assess the likelihood of finding spurious sources in each map. 

As an initial test of the data, each group focused on 
the SXDF data, extracting a list of blank-field, 3 a source 
candidates. From these source lists, a preliminary cross- 
identification candidate source list was created using a simi- 
lar grading scheme to that described in Section ll^ To deter- 
mine the likely number of spurious objects in this list of 'pos- 
itive' sources, each group then applied the same source ex- 
traction methods to inverted, or 'negative' versions of these 
same SXDF maps (cf. Sectionim. The result was that in the 



combined SXDF lists more 'negative' than 'positive' source 
candidates were identified (142 versus 131). From these lists, 
we then removed those source candidates found in noisier re- 
gions of the maps (see Figs. 1 and 2), as these should be less 
reliable. This resulted in an equal number of 'positive' and 
'negative' source candidates in the combined SXDF list (68). 
Higher S/N threshold cuts (up to 5a) did not yield an ex- 
cess of 'positive' over 'negative' source candidates. Spatial 
and temporal tests (see Appendix IB2t were performed 
on Reduction D's > 3a source candidates, and we found 
that the majority of the sources fit within the allowed 2 a 
area of the map's distribution and therefore could not 
be rejected on these grounds. Based on these analyses, no 
single 450 /im blank-field source candidate can be claimed 
as a reliable detection in the SHADES data. 

Although Reductions B and C do agree on the detec- 
tion of a few 450 fim counterparts to the 850 fim sources in 
the LH field (some of which was observed under excellent 
weather conditions as part of the SCUBA 8-mJy survey), in 
no case do all three reductions agree on a detection with a 
consistent position and significance level. As such, we have 
chosen to adopt the 3 a limits on the 450 /im fiux densities 
for subsequent analyses, such as the photometric redshift es- 
timates (Aretxaga et al. in preparation), which benefit from 
the shorter wavelength data. 

B3.3 450 fim photometry 

As the reductions show broad agreement in the 450 ^m fiux 
densities within the uncertainties, and there are no large 
systematic differences, we adopt the 450 fim flux densities 
and photometric errors of Reduction B. We do not claim 
any individual 450 /^m detections. Equivalent 3 a upper lim- 
its are calculated for each source in the following way. We 
construct an error distribution from a histogram of pixel 
S/N. This error function is nearly Gaussian but has slightly 
larger wings. For each source, the error distribution is scaled 
to match the peak and a of the source. We quote in Table|21 
as upper limits the flux densities bounding 99.73 per cent of 
the area of the error function (corresponding to the percent- 
age area between the tails of the 3 a region under a Gaussian 
distribution) . 

B3.4 450 /im stacking analyses 

Although we are unable to claim a 450 fim detection of any 
850 /im SHADES source, we attempt to determine if the 
population as a whole is detected in these data. To do this, 
we perform a series of stacking analyses on the 450 /im data 
at the positions of the 850 fim sources, and also their pro- 
posed 1.4 GHz radio counterparts (Ivison et al. in prepara- 
tion). These analyses were performed independently by each 
group following 2 strategies: (1) pointed 450 /im photometry 
at the precise 850 /tm-selected source positions (this result 
will be biased low since the peak of emission may be offset 
by an amount as great as the JCMT pointing error of 2-3 
arcsec); and (2) a search for the nearest 450 /im peak within 
a 7 arcsec radius of the 850 fim source position, if one exists 
(an estimate which should be biased high) . A 7 arcsec search 
radius was chosen since the probability of a spurious 450 /im 
source lying within this distance of a known 850 /im is very 
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Table Bl. 450 /jm stacked flux densities (in mjy) for the B, C and D reductions. The results are for average flux densities measured at 
the indicated positions, related to those of the 850 /im SHADES catalogue sources. Reduction B has calculated unweighted mean stacked 
flux densities, while Reductions C and D have calculated mean stacked flux densities with inverse noise variance weighting. 



Reduction 


850 fim position 


highest S/N detection 


850 ^m position 


850 ^m position 


radio position 






within Tarcsoc radius of 










(all) 


850 fj,m position (all) 


(radio ID subset) 


(non-radio ID subset) 


(radio ID subset) 


LH 


B 


10.6 ±3.8 


38.3 ±4.0 


16.9 ±4.3 


10.6 ±4.8 


19.2 ±4.2 


C 


12.4 ± 2.1 


30.6 ±2.1 


20.8 ±3.3 


7.4 ± 2.7 


22.0 ±3.3 


D 


10.5 ± 2.3 


28.8 ±2.9 


19.6 ±3.8 


4.9 ± 2.9 


18.8 ±3.8 


SXDF 


B 


5.7 ± 4.6 


24.3 ± 3.7 


11.1 ± 4.2 


14.6 ±5.4 


16.6 ±4.2 


C 


12.5 ±2.4 


37.5 ± 2.4 


14.1 ± 3.9 


11.4 ±3.1 


17.9 ±3.8 


D 


9.9 ±2.7 


36.7 ±3.4 


12.1 ± 4.3 


7.9 ±3.5 


14.0 ±4.2 



low (see lFox et alll2002h . The results from each reduction 
are given in columns 2 and 3 of Table IHTI 

To determine how frequently the stacked flux density 
measurements would occur given the number of 850 fxm 
sources in each field, 10,000 Monte Carlo simulations were 
performed with uniformly-selected random positions in the 
Reduction D maps. A stacked flux density greater than or 
equal to the measured value at the 850 /xm-selected posi- 
tions (strategy 1) occurs < 0.1 percent of the time at ran- 
dom in both fields. The simulations therefore indicate a 
significant detection of the SHADES catalogue of 850 fim- 
selected galaxies at 450 /im. Another indication that we have 
marginally detected the 850 /im sources at 450 fim, is shown 
by the slight positive skewness in the distribution of S/N val- 
ues at the 850 /xm SHADES catalogue positions (as shown 
in Fig. IBll . The simulations were then repeated for strategy 
2, searching for the nearest 450 /im peak within 7arcsec of 
the SHADES catalogue 850 /^m positions. We found that a 
value above the measured stacked value occurred 90 (22) per 
cent of the time in the LH (SXDF) field. The stacked 450 ^m 
flux density obtained using strategy 2 is therefore not high 
enough compared with the Monte Carlos (which pick up 
450 pm noise peaks within the 7 arcsec search radius) to be 
regarded as statistically significant. 

The ratio of stacked 450 fim fiux density to 850 fim flux 
density of SHADES catalogue sources is around 2-2.5, de- 
pending on the precise choice of data and reduction (and 
including the effect of 850 ^m flux deboosting). This ratio 
is low compared to what one would expect, independent of 
any reasonable choice for a redshifted sample of local SED 
templates (which typically yield S'45o/<S'850 <; 4). We found 
that fi tting a range of template gala^xies from lAretxaga et alJ 
l|2003l) yielded redshifts ^ 3 for the population, which is 
higher than that found for SMGs. This result of low stacked 
450 nm fiux density is co nsistent with that found for SMGs 
in the GOODS-N field llPope et al.1 12005') and suggests a 
systematic bias when pushing SCUBA data to these faint 
stacked values. 

The 450 ^m stacking analysis was also carried out at the 
positions of the preliminary sub-sample of radio-identified 
sources (see Section 14.711 . and compared with the results 
of an analysis at the positions of the non-radio-identified 



z. 




450|im S/N Measured at 850nm SHADES Positions 

Figure Bl. Histogram of stacked flux densities in the 450 ^m 
S/N map at the positions of the 850 /^m-selected SHADES cata- 
logue positions (see Table 1^ . The stacked flux densities for the 
LH and SXDF are shown by the large-binned thick histograms. 
The finer-binned histograms represent the 450 fim pixel values 
from the S/N LH and SXDF maps. Notice the appearance of ex- 
cess positive 450 /xm S/N when stacked at the 850 fim positions, 
indicating that we do detect an overlap between these popula- 
tions. While any single 850 /xm object is not reliably detected, 
the 850 fim population as a whole is detected statistically. This 
plot is specifically for Reduction D; the results based on the other 
two reductions are similar. 



sub-sample. Note that the radio-de tected SMGs are mor e 
likely to be detected at 450 /im (e.g. IChapman et alJl20o'5l) . 
The stacked 450 fim flux densities at the 850 fim positions 
of the sources with radio IDs are similar to the stacked 
450 ^m flux densities at the radio positions. This implies 
that moving to the positions of the proposed radio iden- 
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tifications (which one might expect would provide a more 
accurate measure of the true source positions) does not ac- 
tually raise the average 450 /xm flux density significantly for 
the same subset of sources. However, a marked difference was 
found between the stacked 450 /im flux densities of the radio- 
identified subset and the non-radio-identified subset. This 
result may be explained by one of the following: (f ) the non- 
radio-identified subset lies at a higher redshift on average; 
or (2) there are spurious sources in the non-radio-identified 
subset J diluting the 450 ^tm stacked flux d ensity measure- 
ment (llvison et alJ 120021 ICreve et al.ll2004l) . Since there is 
no evidence that the radio-identified SHADES sources are 
more likely to be real, the radio-identified subset is probably 
biased to lower redshifts dChapman et alJl2005^ . yielding a 
higher S'45o/<5'850 ratio. 



APPENDIX C: NOTES ON INDIVIDUAL 
SOURCES 

Here we give notes for some of the individual SHADES 
sources. 

LOCK850.1: This source has a relatively poor spatial 
(3.f3, where the +2 a range of the distribution extends 
to 2.i6), which could perhaps be explained by its proximity 
to LOCK850.4f (22arcsec away). 

LOCK850.2: More than 50 per cent of data in this re- 
gion of the map were afi'ected by the noise spike. Treating 
the noise spike as discussed in Appendix ^resulted in a 2 
per cent decrease in the flux and S/N and no discernible 
positional offset. 

LOCK850.5: Treating the noise spike as discussed in 
Appendix^resulted in an apparent offset in RA of 3 arcsec, 
an offset of 3 arcsec in Dec, and approximate differences in 
flux, noise and S/N of 2, 4, and 2 per cent, respectively. 

LOCK850.8: This source has a relatively poor temporal 
(— i.72, where the — 2(t range of the distribution extends 
to -1.60). 

LOCK850.10: This source was detected in the noisier 
edge region of the map (cr > 3mJy). 

LOCK850.il: Treating the noise spike as discussed in 
Appendix 1X1 resulted in an apparent offset of —3 arcsec in 
Dec. and approximate differences in flux, noise and S/N of 
1, 2, and per cent, respectively. 

LOCK850.12: This source has a relatively poor tem- 
poral x^ (~l-70, where the —2 a range of the distribution 
extends to —1.60). 

LOCK850.U: LOCK850.18 is a nearby neighbour, 
20 arcsec away. 

LOCK850.15: This source has a relatively poor tem- 
poral (~2.12, where the 2 a of the distribution extends 
to —1.60). This source also has a relatively poor spatial x^ 
(7.82, where the -\-2a range of distribution extends to 2.16). 
This source was detected in the noisier edge region of the 
map (a > 4mJy). 

LOCK850.18: LOCK850.14 is a near neighbour, 
20 arcsec away. 

LOCK850.21: Treating the noise spike as discussed in 
Appendix 1^ resulted in an apparent offset of 2.9 arcsec in 
RA and —3 arcsec in Dec, and approximate differences in 
flux, noise and S/N of 1, 1, and 2 per cent, respectively. 

LOCK850.22: Treating the noise spike as discussed in 



Appendix 1X1 resulted in an apparent offset of —3 arcsec in 
RA and approximate differences in flux, noise and S/N of 5, 
2, and 2 per cent, respectively. 

LOCK850.33: Treating the noise spike as discussed in 
Appendix^resulted in an apparent offset of 3 arcsec in RA 
and —3 arcsec in Dec, and approximate differences in flux, 
noise and S/N of 3, 2, and per cent, respectively. 

LOCK850.34'- This source was detected in the noisier 
edge region of the map (cr > 3mJy). 

LOCK850.35: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix^resulted in a 4 per 
cent decrease in the flux, 2 per cent decrease in S/N and no 
noticeable positional offset. This source has a relatively poor 
temporal (—1.94, where the —2 cr range of the distribution 
extends to —1.60). 

LOCK850.39: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^ resulted in a 7 
per cent decrease in the flux, 5 per cent decrease in S /N and 
no noticeable positional offset. 

LOCK850.40: This source has a relatively poor tem- 
poral x^ (—2.05, where the —2a range of the distribution 
extends to —1.60). 

LOCK850.41'- This source has a relatively poor spa- 
tial (5.98, where the +2 a range of the distribution ex- 
tends to 2.16), which could be explained by its proximity to 
LOCK850.1 (22 arcsec away). 

LOCK850.43: Treating the noise spike as discussed in 
Appendix IXI resulted in an apparent offset of —3 arcsec in 
Dec. and approximate differences in flux, noise and S/N of 
12, 0, and 12 per cent, respectively. 

LOCK850.47: Treating the noise spike as discussed in 
Appendix IXI resulted in an apparent offset of —3 arcsec in 
RA and approximate differences in flux, noise and S/N of 6, 
2, and 3 per cent, respectively. 

LOCK850.48: Treating the noise spike as discussed in 
Appendix IXI resulted in an apparent offset of —3 arcsec in 
Dec. and approximate differences in flux, noise and S/N of 
4, 3, and 3 per cent, respectively. 

LOCK850.60: Treating the noise spike as discussed in 
Appendix^resulted in an apparent offset of 3 arcsec in Dec. 
and approximate differences in flux, noise and S/N of 6, 2, 
and 3 per cent, respectively. 

LOCK850.66: Treating the noise spike as discussed in 
Appendix^resulted in an apparent offset of 3 arcsec in Dec. 
and approximate differences in flux, noise and S/N of 2, 1, 
and 5 per cent, respectively. 

LOCK850.67: Treating the noise spike as discussed in 
Appendix^resulted in an apparent offset of 3 arcsec in Dec. 
and approximate differences in flux, noise and S/N of 3, 1, 
and 3 per cent, respectively. 

LOCK850.76: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^resulted in a 13 
per cent decrease in the flux, a 9 per cent decrease in the 
S/N and no discernible positional offset. 

LOCK850.77: Treating the noise spike as discussed in 
Appendix IXI resulted in an apparent offset of —3 arcsec in 
RA and approximate differences in flux, noise and S/N of 1, 
0, and 3 per cent, respectively. 

LOCK850.79: Treating the noise spike as discussed in 



© 0000 RAS, MNRAS 000, 000-000 



36 Coppin et al. 



Appendix 1X1 resulted in an apparent offset of — Sarcsec in 
RA, 3arcsec in Dec, and approximate differences in flux, 
noise and S/N of 4, 4, and 9 per cent, respectively. 

LOCK850. 100: This source was detected in the noisier 
edge region of the map (cr > 3mJy). 

SXDF850. 1 : This source appears to be extended in the 
map in the NS direction. 

SXDF850.5: More than 50 per cent of data in this re- 
gion of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^resulted in a 6 
per cent decrease in the flux, a 2 per cent decrease in S/N 
and no discernible positional offset. 

SXDF850.6: This source has a relatively poor spatial 
(—2.28, where the —2 a range of the distribution extends 
to -2.05). 

SXDF850. 7: This source has a relatively poor spatial 
(2.03, where the +2 a range of the distribution extends to 
2.05). The poor fit could be a result of it being separated by 
less than 22 arcsec from a lower significance rejected source 
(SXDF850.110). 

SXDF850.9: More than 50 per cent of data in this re- 
gion of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^resulted in a 1 
per cent decrease in the flux, and no discernible change in 
S/N or position. 

SXDF850.10: This source has a relatively poor tempo- 
ral (2.27, where the -1-2 a range of the distribution extends 
to 1.60). 

SXDF850.18: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^resulted in a 6 
per cent decrease in the flux, a per cent change in S/N and 
a positional offset of —3 arcsec in Dec. 

SXDF850.19: Treating the noise spike as discussed in 
Appendix^resulted in an apparent offset of 3 arcsec in RA 
and approximate differences in flux, noise and S/N of 1, 2, 
and per cent, respectively. This source has a relatively poor 
spatial (2.75, where the +2 a range of the distribution 
extends to 2.05). 

SXDF850. 21 : This source has a relatively poor spatial 
X^ (2.23, where the -\-2a range of the distribution extends 
to 2.05). This source also has a relatively poor temporal x^ 
(1.84, where the -1-2 cr range of the distribution extends to 
1.60). 

SXDF850.22: Treating the noise spike as discussed in 
Appendix^resulted in an apparent offset of 3 arcsec in RA, 
—3 arcsec in Dec, and approximate differences in flux, noise 
and S/N of 11, 2, and 7 per cent, respectively. 

SXDF850.27: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^ resulted in a 
5 per cent decrease in the flux, a 3 per cent decrease in 
S/N and no noticeable positional offset. This source has a 
relatively poor spatial x^ (2.15, where the +2a range of the 
distribution extends to 2.05). 

SXDF850.28: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix^ resulted in a 10 
per cent increase in the flux, a 13 per cent increase in S/N 
and a positional offset of 3 arcsec in RA. This source has a 
relatively poor temporal x^ (~l-75, where the —2 a range of 



the distribution extends to —1.60). This source was detected 
in the noisier edge region of the map {a > 3 mjy) . 

SXDF850.37: Treating the noise spike as discussed in 
Appendix IXI resulted in an apparent offset of —3 arcsec in 
RA, and approximate differences in flux, noise and S/N of 
1, 1, and per cent, respectively. 

SXDF850.45: This source has a relatively poor tem- 
poral x^ (—1.85, where the —2a range of the distribution 
extends to —1.60). This source was detected in the noisier 
edge region of the map {a > 4mJy). 

SXDF850.47: This source has a relatively poor tempo- 
ral x^ (1.80, where the +2 a range of the distribution extends 
to 1.60). 

SXDF850.50: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^resulted in a 10 
per cent decrease in the flux, a 10 per cent decrease in S/N 
and no discernible positional offset. 

SXDF850.56: Treating the noise spike as discussed in 
AppendixIXIresulted in an apparent offset of 3 arcsec in Dec, 
and approximate differences in flux, noise and S/N of 4, 1, 
and 3 per cent, respectively. 

SXDF850.63: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^ resulted in a 5 
per cent decrease in the flux, no discernible change in S/N 
and a positional offset of —3 arcsec in RA. 

SXDF850.69: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^resulted in a 15 
per cent increase in the flux, a 15 per cent increase in S/N 
and a positional offset of —3 arcsec in RA. 

SXDF850.76: Treating the noise spike as discussed in 
Appendix IXI resulted in an apparent offset of —3 arcsec in 
RA, —3 arcsec in Dec, and approximate differences in flux, 
noise and S/N of 1, 3, and 3 per cent, respectively. 

SXDF850.86: This source has a relatively poor tempo- 
ral x^ (2.21, where the -f 2 a range of the distribution extends 
to 1.60). 

SXDF850.93: Treating the noise spike as discussed in 
Appendix^resulted in an apparent offset of 3 arcsec in Dec. 
and approximate differences in flux, noise and S/N of 3, 1, 
and per cent, respectively. 

SXDF850.95: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^ resulted in a 8 
per cent decrease in the flux, a 6 per cent decrease in S/N 
and positional offsets of —3 arcsec in RA and —3 arcsec in 
Dec. 

SXDF850.96: More than 50 per cent of data in this 
region of the map were affected by the noise spike. Treating 
the noise spike as discussed in Appendix ^ resulted in a 4 
per cent decrease in the flux and a 3 per cent decrease in 
S/N, with no noticeable positional offset. 
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